Methods and systems for particle based treatment using microdosimetry techniques

ABSTRACT

Methods and systems are described for determining particle treatment information. An example method may comprise determining a segment-averaged dose-averaged restricted linear energy transfer. The linear energy transfer may be determined by accounting for variations in segment length of paths of particles in a site.

BACKGROUND

In radiation therapy, the absorbed dose by cells determines theprobability of a cell's survival. However, the biological effectivenessof the same physical dose of particle therapy and x-rays is different¹.In particular, in vivo experiments have shown that protons are notequally effective in the entrance region, the Bragg peak or the distalfalloff region of the beam². This difference in radiobiologicaleffectiveness has been related to the increase of linear energy transfer(LET) towards the end of the proton range³⁻⁵. There is, therefore, anincreasing interest in using LET in the optimization process of thetreatment plan^(6,7). Thus, there is a need for more sophisticatedtechniques for optimizing treatment plans.

SUMMARY

Methods and systems are disclosed for determining treatment information.An example method can comprise determining data indicative of a particlebeam to apply to an object of interest and determining imaging dataassociated with the object of interest. The imaging data can comprise aplurality of voxels of data. The method can comprise determining, foreach of a plurality of domains in a first voxel of the plurality ofvoxels, a plurality of distributions for characterizing interactions ofparticles of the particle beam with the corresponding domain. Theplurality of distributions can comprise a first distribution indicativeof energy imparted in a corresponding domain due to a particle travelingin the domain. The method can comprise determining, based on theplurality of distributions and for each of the plurality of domains,first data comprising energy imparted by the particle beam to thecorresponding domain. The method can comprise outputting, based on thefirst data, data associated with treatment of the object of interest bythe particle beam.

An example method can comprise determining volumetric data associatedwith an object of interest. The volumetric data comprise a plurality ofdomains. The method can comprise determining, for at least a portion theplurality of domains, a plurality of distributions for characterizinginteractions of particles of a particle beam with the correspondingdomain. The plurality of distributions comprises a first distributionindicative of energy imparted in a corresponding domain due to aparticle traveling in the domain. The method can comprise determining ananalytical function based on one or more of the plurality ofdistributions. The method can comprise determining, based on theanalytical function and for each of the plurality of domains, first datacomprising energy imparted by the particle beam to the correspondingdomain. The method can comprise outputting, based on the first data,data associated with treatment of the object of interest by the particlebeam.

This Summary is provided to introduce a selection of concepts in asimplified form that are further described below in the DetailedDescription. This Summary is not intended to identify key features oressential features of the claimed subject matter, nor is it intended tobe used to limit the scope of the claimed subject matter. Furthermore,the claimed subject matter is not limited to limitations that solve anyor all disadvantages noted in any part of this disclosure.

BRIEF DESCRIPTION OF THE DRAWINGS

The accompanying drawings, which are incorporated in and constitute apart of this specification, illustrate embodiments and together with thedescription, serve to explain the principles of the methods and systems.

FIG. 1A shows an overview for calculating treatment information.

FIG. 1B shows different segment lengths for particles in a domain.

FIG. 1C shows an overview for calculating treatment information.

FIG. 1D shows a general technique for calculating treatment information.

FIG. 2A shows a first part of an example technique for calculatingtreatment information.

FIG. 2B shows a second part of the example technique for calculatingtreatment information.

FIG. 3 shows a setup for simulations carried out in Geant4-DNA.

FIG. 4A shows graphs of calculations from MC simulations andcorresponding fitted function for average single-event energy impartedto sites with diameters of 1 μm, 5 μm, and 10 μm.

FIG. 4B shows graphs of calculations from MC simulations andcorresponding fitted function for the standard deviation of thecorresponding distributions of energy imparted to sites with diametersof 1 μm, 5 μm, and 10 μm.

FIG. 5 shows graphs of calculations from MC simulations andcorresponding fitted function for average segment length in sphericalsites with diameters of 1 μm, 5 μm, and 10 μm.

FIG. 6 shows graphs of dose-mean lineal energy obtained from MCsimulations and with proposed analytical models (lines) for sphericalsites with diameters of 1 μm, 5 μm, and 10 μm.

FIG. 7 shows graphs comparing the analytical expression shown inequation (13) and the MC data resulting from Geant4-DNA simulations forthe average and the standard deviation of the distribution of impartedenergy per collision.

FIG. 8A show shows graphs of segment-averaged dose-averaged restrictedLET for spherical site with diameters of 1 μm, 5 μm, and 10 μm.

FIG. 8B shows absolute differences between segment-averageddose-averaged restricted LET calculated with the proposed analyticalmodels and with Geant4-DNA MC simulations.

FIG. 9A shows frequency-average of single-event energy imparted to thesite for sites of 3 μm and 7 μm, interpolating from the data for sitesof 1 μm, 5 μm and 10 μm.

FIG. 9B shows segment-averaged dose-averaged restricted LET for sites of3 m and 7 μm, interpolating from the data for sites of 1 μm, 5 μm and 10am.

FIG. 10 is a block diagram illustrating an example computing device.

FIG. 11A shows residuals between MC-based and model-based calculationsof average single-event energy imparted to the site.

FIG. 11B shows residuals between MC-based and model-based calculationsof the standard deviation of the corresponding distributions of energyimparted

FIG. 12A shows residuals between MC-based and model-based calculationsof the average segment length for spherical sites with diameters of 1am, 5 μm and 10 am, respectively.

FIG. 12B shows residuals between MC-based and model-based calculationsof the weighted-average segment length for spherical sites withdiameters of 1 am, 5 μm and 10 am, respectively

FIG. 13 shows calculations from MC simulations done with Geant4-DNA inliquid water for monoenergetic protons from 0.01 to 100 MeV (points) andcorresponding fitted function, equation (12) (lines) for theweighted-average segment length described in spherical sites withdiameters of 1 μm, 5 μm and 10 μm, respectively

FIG. 14 shows residuals between the dose-mean lineal energy obtainedfrom MC simulations of monoenergetic protons in liquid water done withGeant4-DNA and with our proposed analytical models by using equation(2), for spherical sites with diameters of 1 μm, 5 m and 10 μm,respectively

FIG. 15A shows calculations for the average of the distribution ofenergy imparted per collision in spherical sites with diameters of 1 μm,5 μm and 10 μm, respectively.

FIG. 15B shows the standard deviation of the distribution of energyimparted per collision in spherical sites with diameters of 1 μm, 5 μmand 10 μm, respectively

FIG. 16A shows residuals between MC-based and model-based calculationsof average of the distribution of energy imparted per collision inspherical sites with diameters of 1 μm, 5 μm and 10 μm, respectively.

FIG. 16B shows residuals between MC-based and model-based calculationsof standard deviation of the distribution of energy imparted percollision in spherical sites with diameters of 1 μm, 5 μm and 10 μm,respectively.

FIG. 17 shows residuals between the weighted average of the distributionof energy imparted per collision, δ₂, calculated for monoenergeticprotons in liquid water from Geant4-DNA MC simulations and from theanalytical models described in equation (13) for spherical sites withdiameters of 1 μm, 5 μm and 10 μm, respectively

FIG. 18A shows residuals between MC-based and model-based calculationsof frequency-average of single-event energy imparted to the site, E, forsites of 3 μm and 7 μm interpolating from data for sites of 1 μm, 5 μmand 10 μm.

FIG. 18B shows residuals between MC-based and model-based calculationsof segment-averaged dose-averaged restricted LET L _(Δ,D,s), for sitesof 3 μm and 7 μm interpolating from data for sites of 1 μm, 5 μm and 10μm. Statistical uncertainties propagated from the MC simulations (1σ)are shown with error bars.

FIG. 19A shows radionuclides uniformly distributed throughout themembrane of a spherical cell.

FIG. 19B shows the distance for each particle coming out of a membranepoint is the distance to a plane perpendicular to the particle track andtangential to the nucleus.

FIG. 20 shows a spectrum ϕ_(E)(E) of alpha particles coming to thespherical target shown in FIG. 19B.

FIG. 21A shows example analytical functions for mean energy imparted tospherical sites with diameters of 1 μm (top), 5 μm (center) and 10 am(bottom), respectively.

FIG. 21B shows example analytic functions for variance of the energyimparted to spherical sites with diameters of 1 μm (top), 5 μm (center)and 10 μm (bottom), respectively.

FIG. 22A shows analytical calculations and Geant4-DNA simulations for aspherical cell of radius 7.5 μm and variable nucleus radius for y _(F)and y _(D).

FIG. 22B shows analytical calculations and Geant4-DNA simulations for aspherical cell of radius 7.5 μm and variable nucleus radius for z _(1,F)and z _(1,D).

DETAILED DESCRIPTION OF ILLUSTRATIVE EMBODIMENTS

For at least the reasons explained the background section, there is needfor an accurate and fast method to calculate LET. The present workproposes the use of generalized microdosimetry theory to calculate LETfor particle beams based on the distribution of energy imparted inmicroscopic (and possibly biologically relevant) structures⁸⁻¹⁰.

The present disclosure introduces the concept of segment-averaged linearenergy transfer (LET) as a new approach to average distributions of LETof particle beams (e.g., proton beams) based on a revisiting ofmicrodosimetry theory. The concept of segment-averaged LET is then usedto generate an analytical model from Monte Carlo simulations data toperform fast and accurate calculations of LET distributions for particlebeams.

The distribution of energy imparted by a particle beam into arepresentative biological structure or site is influenced by thedistributions of (1) LET, (2) segment length, which is the section ofthe proton track in the site, and (3) energy straggling of the particlebeam. The distribution of LET is thus generated by the LET of eachcomponent of the beam in the site. However, the situation when the LETof each single proton varies appreciably along its path in the site isnot defined. Therefore, a new distribution can be obtained if theparticle track segment is decomposed into smaller portions in which LETis roughly constant. The term “segment distribution” of LET can be theone generated by the contribution of each portion. The average of thatdistribution is called segment-averaged LET. This quantity is obtainedin the microdosimetry theory from the average and standard deviation ofthe distributions of energy imparted to the site, segment length andenergy imparted per collision. All this information is calculated forprotons of clinically relevant energies by means of Geant4-DNAmicrodosimetric simulations. Finally, a set of analytical functions isproposed for each one of the previous quantities. The presented modelfunctions are fitted to data from Geant4-DNA simulations formonoenergetic beams from 100 keV to 100 MeV and for spherical sites of 1μm, 5 μm and 10 μm in diameter. The same approach can be employed forsites of different shapes and sizes.

The results described further herein show that the average differencesalong the considered energy range between calculations based on thedisclosed analytical models and MC for segment-averaged dose-averagedrestricted LET are −0.2±0.7 keV/μm for the 1 μm case, 0.0±0.9 keV/μm forthe 5 μm case and −0.3±1.1 keV/μm for the 10 μm case, respectively. Allaverage differences are below the average standard deviation (1σ) of theMC calculations.

The presently disclosed techniques comprise a novel approach foraveraging LET for a particle beam to incorporate the effects produced bythe variation of stopping power of each individual proton alongmicroscopic biological structures. An analytical model based on MCsimulations allows for fast and accurate calculations ofsegment-averaged dose-averaged restricted LET for particle beams, whichotherwise would need to be calculated from exhaustive MC simulations ofclinical plans.

The techniques and models disclosed herein can be used to calculate dose(e.g., simultaneously with LET) by using specific energy (z) instead oflineal energy (y). This calculation of dose can be performed by dividingthe energy imparted by the mass of the site (e.g., instead of thesegment length distribution). Using the techniques and/or modelsdisclosed herein both averages and weighted averages for z and y can bedetermined. Any biologically relevant quantity derived from dose, LET,specific energy, and/or lineal energy can be determined (e.g., andoutput as part of a treatment planning system).

FIGS. 1A-D provide an overview of the disclosed methods and systems.FIG. 1A shows an overview for calculating treatment information. Aparticle beam can strike a biological structure (e.g., in a site,domain). The present technique adds a new approach to usingmicrodosimetry to determine the energy imparted to a microscopicstructure. Unlike conventional techniques, the present approach accountsfor particles (e.g., protons) stopping in a site and for scenarios whenLET changes in a site. The actual segment length of the particle pathcan be accounted for in the present technique. The microscopic quantitylineal energy y may be determined using the equation shown, which isdescribed further herein. The lineal energy may be used to calculate themacroscopic quantity of LET.

FIG. 1B shows different segment lengths for particles in a domain. Anobject of interest may be scanned and/or imaged to produce imaging data(e.g., or scanning data). The imaging data may be two-dimensional orthree-dimensional. The imaging data may comprise a plurality of voxels(e.g., each voxel as a single unit of data). A single voxel may besubdivided into a plurality of domains (e.g., or sites). A domain (e.g.,or site) may be a geometrically defined volume in a voxel. FIG. 1B showsan example domain. A domain (e.g., or site) may have one or morestructures of the object of interest. The domain may be modeled as asphere as shown, but other shapes may also be used. As shown, dependingon the particle, the particle may have a different segment length. Thesegment length may comprise a distance a particle travels from enteringthe domain until coming to rest (e.g., or exiting the domain).

FIG. 1C shows an overview for calculating treatment information.Restricted LET can be dependent on the energy that secondary electronsneed to escape the domain. Site-average LET can depend on the length ofthe proton track in the site. This figure shows an example diagram forMonte Carlo simulations with G4-DNA. The particle track (e.g., shown asa line into the sphere) may be simulated starting from the same fixedpoint (e.g., a source). Then, the site position may be sampled for eachevent at a certain point in the gray region shown, so that the relativeposition between the particle track and the sphere varies with theevent.

FIG. 1D shows a general technique for calculating treatment information.As shown, LET can be calculated for a site/domain of dimension given bydiameter d. A general equation is shown for calculating LET for anybeam. Also shown, is that for the particular case of LET not changingand particles not coming to rest in a domain, Kellerer's equation may bederived.

FIGS. 2A-B shows an example technique for calculating treatmentinformation. FIG. 2A shows a first part of an example technique forcalculating treatment information. FIG. 2B shows a second part of theexample technique for calculating treatment information. This disclosedprocess may be used to determine LET, dose, and/or other treatmentrelated information. The top box shows a sample calculation for a domain(e.g., or site) of dimension d1.

Data indicative of a particle beam may be determined. The dataindicative of the particle beam many comprise data from a treatmentplanning system, data indicative of a particle beam to be generated,and/or the like. The characteristics of the particle beam coming fromthe nozzle are usually given by the treatment planning system (TPS).This software comprises models for the beam that are fitted toexperimental data obtained for each single treatment machine/room.Therefore, from the TPS can be obtain: (I) the number of protons comingto the patient surface, (II) the angular aperture of the beam, (III) theenergy spectrum, a combination thereof, and/or the like.

The data indicative of the particle beam may comprise data specific toeach voxel of patient data. The data indicative of the particle beam maycomprise data indicative of a number of protons, an angular aperture ofbeam, and/or an energy spectrum may be determined at each voxel of a 3Dpatient data set.

The data indicative of the particle beam may be determined by using asimulation (e.g., or modeling) process, such as a Monte Carlo code(e.g., much faster and simpler than the microdosimetric). The simulationprocess can be configured to determine (e.g., for each of a plurality ofdomains, sites, voxels, points, and/or the like) how many protons arepredicted to be delivered from a beam, how broad is predicted beam,and/or the spectrum of the beam.

A plurality of distributions may be determined. The plurality ofdistributions may be for characterizing interactions of particles of aparticle beam with a corresponding domain (e.g., or site). The pluralityof distributions can comprise a first distribution indicative of segmentlength of a particle path in the site, a second distribution indicativeof an energy imparted in the domain due to a particle coming to rest,and a third distribution indicative of energy imparted in the site dueto collisions of a particle. Determining the plurality of distributionscan comprise using one or more simulations to generate the plurality ofdistributions. Example simulations may comprise Monte Carlo simulations,GEANT4-DNA simulations, and/or the like. GEANT4-DNA can be a particularMonte Carlo code. The simulation can comprise using random numbers todetermine a value from a given range of possibilities. The simulationcan consider the processes possibly occurring in the interactionparticle-matter (whose relative probabilities are known) and selectingone of possibility by generating a random number. By simulating asufficiently large number of cases (e.g., particles in a beam) theresulting distribution can provides an estimation of a real-worldresult. The simulation may simulate transport of particles in matter.The simulation can simulate interactions between particles and matter ata microscopic level.

An energetic kernel can be determined based on the plurality ofdistributions. The energetic kernel can comprise a first average of atleast one of the plurality of distributions and a first variance of theat least one of the plurality of distributions. For example, an averageof the distribution and a variance of the distribution can be determinedfor each of the plurality of distributions. This determination can bemade for each proton energy of an expected particle beam.

The term “energetic kernel” indicates that the model can transform agiven function of the particle energy (e.g., the energy spectrum) intoanother function independent of the particle energy (e.g., dose or LET)by integrating over the particle energy the product of the input and theenergetic kernel (e.g., which can also be a function dependent on theparticle energy). Continuing at FIG. 2B, the averages and the variancescalculated (e.g., the energetic kernel) may each undergo a convolutionwith the energy fluence (e.g., for a particular particle beam, such as aproton beam). The results of the convolution may be used to calculatethe segment-averaged dose-averaged restricted LET for a given sitesize/diameter. The result can be used to further refine the propertiesof the particle beam.

The results can be used to determine (e.g., or update a treatment plan).The treatment plan can be a plan to deliver a dose to correspondingvoxels of a set of voxels in a 3D image that represent the patient. Apredicted distribution of dose and/or a distribution of LET can bedetermined. The predicted distribution of dose and/or a distribution ofLET can comprise a value of the corresponding quantities for each of thevoxels. The treatment plan can comprise delivery of more than one beam.This is done for several reasons: to produce a distribution of doseuniform and conformed to the tumor, to distribute the dose to healthytissue on different entries (e.g., so that each entry carries a lowerdose but all of them converge in the tumor), or to decrease theuncertainties. Using the techniques described herein, a distribution ofLET coming from each beam can be calculated. A distribution of LET forthe whole plan can also be calculated (e.g., based on the distributionfor each beam). The treatment plan can comprise one or more treatmentparameters. The treatment parameters can be determined, adjusted,updated, modified, optimized, and/or the like based on a determineddose, LET, energy imparted to a domain, any other value calculated usingthe techniques disclosed herein, and/or the like. The treatmentparameters can comprise a number of beams to be used, angles ofincidence, energy and field size for each one of the beams, relativeweight (which is, roughly, the number of particles coming from thatbeam) of the beam, and/or the like.

The results of the techniques disclosed herein can comprise adetermination of values for each of a plurality of sites/domains. Avoxel (e.g., normally, much larger than a site) can be considered ascompletely filled by a regular grid of sites/domains. The number ofprotons in the voxel can be determined. The number of events produced bythose protons in all the sites/domains contained in the voxel can bedetermined. The results can be divided by the number of sites in thevoxel.

Additional details are described using the follow detailed explanationand examples. The concept of electronic stopping power—or unrestrictedLET—could be considered to measure the energy deposited by particlebeams due to electronic collisions in cellular or sub-cellular regions.However, there are at least two relevant considerations disregarded bythis quantity: the energy transported out of these regions by secondaryelectrons and the finite range of protons¹¹; in other words, situationsof fast changes of LET in short paths, which occurs in the trajectoryend of the beam. The first can be incorporated in the concept ofrestricted linear energy transfer⁸ L_(Δ), where the energy deposition bysecondary electrons with kinetic energy beyond a cut-off value Δ isdisregarded. Efforts have been devoted to calculate either unrestrictedor restricted LET for clinical particle beams in acceptable times inclinical routine. Some analytical approaches have been proposed based onthe proton stopping power values¹²⁻¹⁵, while other works try to includethe effect of the energy carried away by secondary electrons usinggeometrical considerations to obtain spatially or radially restrictedLET¹⁶⁻¹⁸. However, how to deal with situations in which LET may changesignificantly across microdosimetric structures remains unclear. This isespecially important since in such situations the highest values of LETare produced.

In this disclosure, calculation method of LET based on microdosimetry isproposed to address simultaneously the limitations referred above, e.g.,both situations in which restricted LET applies and when LET variesalong microscopic structures. The resulting distributions from thesesimulations were used in this work as input data to generate analyticalmodels to calculate LET, and the following sections describe thesemodels in detail. The use of these models to calculatemicrodosimetry-based LET in clinical particle beams may be the onlyaccurate approach available in clinically reasonable times, since fullMonte Carlo microdosimetric LET calculations would need to descend tothe track-structure level of detail, which is unaffordable in terms ofcomputation.

Example methods and materials in accordance with the present disclosureare described as follows.

Definitions of relevant LET and microdosimetry concepts are described asfollows.

Linear energy transfer (LET) is a non-stochastic quantity that considersthe mean energy lost by a charged particle per unit path length. Whendealing with a field of charged particles with diverse energies, as inthe case of a clinical beam, a different value of LET can be found foreach particle. Therefore, rather than a single LET value, a spatialdistribution of LET has to be considered in a patient. Thus, it isconvenient to reduce that distribution to average values, which can bedone in two different ways: one is the fluence-based average,traditionally called track-averaged LET, L _(T); while in the secondone, is called dose-averaged LET, L _(D), where the LET of each particleis weighted by its contribution to the local dose¹². When energyimparted into a certain region is considered instead of the energy lostby the particle, the concept of energy-restricted LET, L_(Δ), becomesuseful.

In microdosimetry, a site is a region of the irradiated volume with agiven shape and microscopic dimension¹⁹. Although any shape can beconsidered, this work is restricted to spherical sites for simplicity.When a charged particle strikes a site, a stochastic amount of energy,E, is imparted to it, carried by both the particle itself and thesecondary electrons generated in electronic collisions. The term“imparted” means that energy is considered to be deposited locally,i.e., within the site. Each energy deposition resulting from a differentprimary particle and its showers of secondary particles is called anevent. In calculating L _(D), only electrons may be taken intoconsideration as secondary particles in each event; this implies thatother types of secondary particles (e.g. alphas) must then be consideredas separate contributors to L _(D). Therefore, when a number ofparticles impact on a site, a distribution of energy imparted per event,so-called single-event distribution²⁰, is originated, with averageenergy imparted ε and variance of σ_(ε) ².

The techniques described herein relates to a clinical particle beamirradiating a biological target. Tissues can be modeled as composed ofsites, also called domains²¹, which may represent cells or sub-cellularcritical structures. The distributions of energy imparted in sites by abeam of particle can be calculated with a given energetic spectrum. Thisbeam is considered to be broad enough to laterally cover all the siteplus a certain margin used to ensure that all energy transfers impartedby secondary electrons generated outside the domain are accounted for.

In this case, for any proton randomly selected from the incidental beam,the portion of its track inside the site is also a stochastic quantity.Protons are assumed to travel in a straight line, which is a goodapproximation as long as non-electronic collisions are disregarded inthe analysis²². The path of the track the proton travels within the siteis called segment length²³ (s). This definition applies to both whenprotons stop in the site and when they completely traverse it. Thisquantity is used to define lineal energy (y_(s)) as the quotient betweenthe energy imparted to the site in an event and the mean segment length:

$\begin{matrix}{y_{s} = {\frac{ɛ}{\overset{\_}{s}}.}} & (1)\end{matrix}$

The subscript s indicates that segment length is used to define linealenergy. This definition of lineal energy is different and novel incomparison to conventional approaches. For example, the terms “segmentlength” and “chord length” can be understood to have different meanings.The segment length can be the actual segment of the particle track inthe site in all cases (e.g., including the cases in which the particlestops in the site, called “stoppers”). Chord length can refer to thegeometrical problem of calculating the distribution of intersectionsbetween a (infinite) straight line and the site (a sphere, in our case).Chord length only takes into consideration particles that completelytraverse a domain or site (called “crossers”). The traditionaldefinition of lineal energy uses chord length instead of segment length,which means that its applicability was restricted only to “crossers”. Byusing segment length, the present techniques extending that definitionto include “stoppers”.

Lineal energy is, therefore, a stochastic quantity with its owndistribution of frequencies, characterized by the frequency probabilitydensity f(y_(s)), whose average is denoted as y _(F,s). A relateddistribution of y_(s) considers the dose absorbed by the tissue in eachevent with lineal energy y_(s) This is called dose distribution ofy_(s), d(y_(s)), and its average is called dose-mean lineal energy²⁰,y_(D,s). Since lineal energy is proportional to the dose absorbed by thesite, the dose distribution of y is essentially a weighted distributionmultiplied by a constant factor, i.e., d(y_(s))∝y_(s)f(y). Therefore,the averages are related by the following equation, which applies to allfrequency and weighted averages⁹:

$\begin{matrix}{{{\overset{\_}{y}}_{D,s} = {{{\overset{\_}{y}}_{F,s}( {1 + \frac{\sigma_{y_{s}}^{2}}{{\overset{\_}{y}}_{F,s}^{2}}} )} = {\frac{\overset{\_}{ɛ}}{\overset{\_}{s}}( {1 + \frac{\sigma_{ɛ}^{2}}{{\overset{\_}{ɛ}}^{2}}} )}}},} & (2)\end{matrix}$

where σ_(y) _(s) ², is the variance of the frequency distribution ofy_(s). The second equality is obtained from the definition of y_(s) inequation (1).

Segment distribution of LET is described as follows.

According to equation (2), y _(D,s) is not only related to the averageof the frequency distribution of lineal energy in the site, but it alsotakes into account the variability of that distribution. It can beshown⁹ that three factors are involved in the different values of y_(s)for each case: the different LET of the protons in the site, thedifferent chord length of each proton and the energy straggling. Thiscan be formulated in terms of the variance of each variability source.By considering equation (2) one can obtain the variance of thedistribution of segment length in terms of the average of its weighteddistribution d(s)=sf(s), also called weighted-mean segment length, s_(D). Additionally, the variance of the straggling distribution can beexpressed as δ₂/ε, where δ₂ is also the weighted-average of thedistribution of energy imparted into the site per collision²⁴.

Disclosed herein is an extension of the applicability of the LET conceptbased on the following main idea: the variability of LET in the site canbe decomposed into two different terms. In the situation in whichprotons have no appreciable change of LET in the site, a field ofprotons striking the site with different LET, or, equivalently,different energy, produce a variety of LET values. Thus, a “fielddistribution” of LET is produced by different protons. On the otherhand, if a single proton loses an appreciable amount of its kineticenergy inside the site, its stopping power and, consequently, its LETare not unique as they can vary considerably in the proton path withinthe site. This case in which a proton changes its L_(Δ) within the sitehas not been previously considered in the microdosimetry literature,since LET is a macroscopic concept and such a change would have astochastic nature. In order to overcome the indetermination on the LETconcept, a new distribution of LET values can be defined based on thisvariability within the site. This distribution can be understood asfollows. Let dx be an element of length small enough so that L_(Δ) canbe considered constant along dx for each proton energy. Then, thedifferent values of L_(Δ) in each dx of the segment length within thesite result in a new LET distribution. This distribution can be called“segment distribution” of LET for each proton traversing the site. This“segment distribution” is also a source of variability in y_(s), sincedepending on the actual track segment that intersects the site, thecombination of L_(Δ) values for the proton in the site might change fromevent to event. Hence, the total distribution of LET in the consideredsegment length is produced by the combination of the “fielddistribution” and the “segment distribution”.

Based on the introduction of a new level of variability in the linealenergy due to the segment distribution defined above, thesegment-distribution-average value of the dose-weighted distribution ofLET is here called “segment-averaged dose-averaged LET”, L _(D,s). This“segment-average” depends on the shape and size of the site. Again, thisquantity is related to the variance of the LET distribution by means ofequation (2). Considering all above, from the relation between relativevariances of the sources of y_(s) variability, the following equationfor the segment-averaged dose-averaged restricted LET can be obtained byemploying analogous arguments to those developed by Kellerer⁹:

$\begin{matrix}{{{\overset{¯}{L}}_{\Delta,D,s} = {{\frac{\overset{\_}{s}}{{\overset{\_}{s}}_{D}}{\overset{¯}{y}}_{D,s}} - \frac{\delta_{2}}{{\overset{\_}{s}}_{D}}}}.} & (3)\end{matrix}$

Here, as energy imparted into a spherical site is considered instead ofenergy lost, the A subscript is added to L _(D,s) in equation (3) toindicate dose-averaged restricted LET. Now, A is equal to the energy forsecondary electrons that makes their range match the distance they haveto travel to escape from the site. Therefore, by employing equation (3),it is possible to obtain a calculation for LET that does not count theenergy escaped from each site or biological domain by secondaryelectrons and the effects of the energy straggling. Additionally, thisdefinition of LET applies to the situation in which the energy impartedper unit length varies appreciably within the site or the extremesituation in which protons stop inside the site.

A particular case of interest: constant LET in the site is described asfollows.

When protons completely traverse the site, i.e., when their range isgreater than the site diameter, the calculation of the distribution ofsegment length becomes a purely geometrical problem and it is obtainedby determining the intersection between a straight line and the site.This is a well-known and solved problem for a spherical site²⁵. In thiscase, the ‘segment length’ is called ‘chord length’ (1). For sphericalsites with diameter d, the chord length distribution is given byf(l)=21/d² (0≤l≤d), with average l=2d/3. This quantity is traditionallyused in microdosimetry to define lineal energy (y) instead of the moregeneral segment length²⁰ defined above. Therefore, in this particularcase, (1) becomes:

$\begin{matrix}{y = {\frac{ɛ}{\overset{\_}{l}}.}} & (4)\end{matrix}$

The case in which protons have a kinetic energy large enough to considerthat their LET does not change within the site is described as follows.In this case, furthermore, protons always cross the site completely thuss→l. Also, under constant LET conditions, only the “field distribution”of LET plays a role. Therefore, in this limit, L _(Δ,D,s) L_(Δ,D) andequation (3) becomes:

$\begin{matrix}{{{\overset{¯}{L}}_{\Delta,D} = {{{\frac{\overset{\_}{l}}{{\overset{\_}{l}}_{D}}{\overset{¯}{y}}_{D}} - \frac{\delta_{2}}{{\overset{\_}{l}}_{D}}} = {{\frac{8}{9}{\overset{\_}{y}}_{D}} - \frac{4\delta_{2}}{3d}}}},} & (5)\end{matrix}$

where the last equality is valid only for spherical sites with diameterd. Equation (5) is not new^(8,9) but now it is presented as a particularcase of Equation (3) to determine the dose-averaged LET of particlebeams. Therefore, its validity can be extended by presenting a moregeneral theory in equation (3), in which this result is included.

Monte Carlo simulations for monoenergetic protons and its application topoly-energetic (clinical) beams are described as follows.

The segment-averaged LET definition is, generally, inconvenient since itdepends on the site shape and dimension and, therefore, producesdifferent values for each site considered. Moreover, segment length forlow energy protons are complicated to calculate analytically, thus usingequation (3) instead of equation (5) seems to add complexity to LETcalculations. However, by using Monte Carlo track structure (MCTS)simulations, all these quantities can be obtained. Furthermore, theproblem of calculating y _(D.s) and L_(Δ,D,s) in a poly-energetic(clinical) beam can be decomposed into simulations of monoenergeticprotons, from which an energy kernel can be developed starting from thecorresponding distributions of energy imparted. Then, if the protonclinical beam is composed by an energy spectrum given by ϕ(E), as theenergy imparted distribution is independent for each monoenergeticproton, the average, ε _(ϕ) and the variance of the total distributionof energy imparted, ϕ_(ε) _(ϕ) ², (i.e. considering the whole fluence)are given, respectively, by

$\begin{matrix}{{{\overset{¯}{ɛ}}_{\phi} = \frac{\int{{\phi(E)}{\overset{¯}{ɛ}(E)}dE}}{\int{{\phi(E)}dE}}},} & (6) \\{{\sigma_{ɛ_{\phi}}^{2} = \frac{\int{{\phi(E)}{\sigma_{ɛ}^{2}(E)}dE}}{\int{{\phi(E)}dE}}},} & (7)\end{matrix}$

where ε(E) and σ_(ε) ²(E) are the average and the variance of thedistribution of energy imparted for a proton of energy E, respectively.The average s and the variance σ_(s) ² of the distribution of segmentlength (which includes the case in which chord length is applicable)might be obtained in a similar way from the monoenergetic cases.Therefore, the y _(s,D) value for a poly-energetic beam can becalculated from equation (2), in which the values calculated usingequations (6-7) have to be inserted for the average and the variance ofthe distribution of energy imparted.

The last quantity to calculate L _(Δ,D,s) by means of equation (3) isδ₂, as s _(D) can be obtained from s and of σ_(s) ². This is theweighted average of the distribution of energy imparted per collision inthe site. Consequently, it can be expressed as:

$\begin{matrix}{{\delta_{2} = {\delta_{1}( {1 + \frac{\sigma_{\delta}^{2}}{\delta_{1}^{2}}} )}},} & (8)\end{matrix}$

where δ₁ and σ_(δ) ² are, respectively, the average and the variance ofthe frequency distribution of energy imparted per collision. These twolatter quantities are additive in the same way that the average and thevariance of the distribution of energy imparted per event are.Therefore, both δ₁ and σ_(δ) ² can be obtained for a poly-energetic beamfrom monoenergetic data by expressions similar to equations (6-7). Ingeneral, L _(Δ,D,s) can be derived for mono-energetic or poly-energetic(clinical) beams using equations (4-5, 8) by collecting monoenergeticdata for ε(E), σ_(ε) ²(E), s(E), s _(D)(E), δ₁(E) and σ_(s) ²(E) andthen integrating over the different energies involved in the spectralfluence distribution of the beam. Additionally, specific energy isdefined as z=ε/m where m is the mass of the site- and its average is theabsorbed dose in the site, D. Therefore, by using equation (6) ispossible to calculate the absorbed dose for the given beam as D=z=ε_(ϕ)/m.

For the reasons indicated in the paragraph above, mono-energetic protonsand their secondary electron tracks were simulated in water usingGeant4-DNA²⁶⁻²⁹. While the description below and elsewhere herein mayreference protons, it should be understood that the same techniques maybe applied to other particles used in a therapy beam. Three differentsite diameters were used to sample microdosimetric energy distributions:1 μm, 5 μm and 10 μm. Proton tracks originated from a point source werearranged to penetrate a water-made box and the position of the sphericalsite was uniformly sampled in a slab of thickness equal to its diameter,which is indeed the length of the track segment analyzed (FIG. 3). Amargin equal to the maximum range of the secondary electrons generatedby protons was added both upstream and downstream the scoring region,where the site is randomly placed, in order to assure charged particleequilibrium condition in the site. This approach is geometricallyequivalent to the situation in which the site is fixed in space and abroad beam traverses it. Geant4-DNA allows the tracking of electronsdown to ˜10 eV, in which case the energy is considered to be depositedlocally. Although protons with energy down to 10 keV have beensimulated, calculations below 100 keV were only taken into account toobtain the distribution of segment length, and not for the distributionof energy imparted. This is because most of the secondary electronsgenerated by protons of such low energies have energies below the lowerlimit of Geant4-DNA, mentioned above. However, these simulations arestill valid for obtaining the distribution of segment length.Furthermore, simulations were performed for proton energies up to 100MeV, which is the upper limit for Geant4-DNA to track protons. To obtainthe distribution of energy imparted per collision, and its solecontribution to the variance of energy imparted, the site position canbe fixed in a way that the proton path intersected it exactly defining alength equal to the mean segment length. Simulations were performed byusing our computing cluster hosted at CICA (Seville, Spain), consistingof 24 computational nodes with 2×12C AMD Abu Dhabi 6344 (2.6 GHz/16 MBL3).

FIG. 3 shows a setup for simulations carried out in Geant4-DNA.Monoenergetic protons emerge from a point source, travelling through awater box. The region in which the site can be sampled (in gray color)is one site diameter in-width and it is laterally extended to ensure anysecondary electron generated by the proton stops in the box and itsenergy can be accounted. Margins upstream and downstream are added toprovide electronic equilibrium conditions in the scoring region. Somepossible site positions are shown in the diagram.

An example model for the average and variance of the distribution ofenergy imparted for monoenergetic proton beams is described as follows.

Analytical functions to model the behavior of the average and thevariance—or, equivalently, the standard deviation—of the distribution ofenergy imparted by protons of different energies are presented in thissection. A phenomenological function of the proton incidental energy forthe average energy imparted (ε(E)) can be defined as follows: it isexpected that the energy imparted is proportional to the stopping powerof the proton at that energy. According to the Bethe formula³⁰ in itsnon-relativistic approximation, this dependency is proportional tolog(a·E)/E, where a=4m_(e)/m_(p)I=27.9 MeV⁻¹, being m_(e) and m_(p) theelectron and proton mass at rest, respectively, and I is the ionizationpotential in water: 78 eV³¹. However, this dependency has a saturatedbehavior due to the finite dimension of the site, e.g., due to thenature of a restricted LET. When the incidental proton has a very lowenergy, all the secondary electrons deposit their energy locally withinthe site and the amount of energy imparted increases as the protonenergy does. However, as energy further increases, some secondaryelectrons can leave the site so that the energy imparted to the site percollision does not increase anymore. This can be modeled by using anerror function. Thus, proposed is ε(E)∝erf(kE^(q))·log(27.9E)/E·R/(R+2d/3), where k and q were introduced as two degrees of freedomto fit the error function to the Geant4-DNA calculations. Furthermore,when proton range is considerably shorter than the site dimension, theprobability for a proton track located laterally far from the sitecenter to enter the site starts to decrease due to the sphere curvature.Therefore, ε(E) should also be proportional to this probability ofinteraction, since the greater the probability, the greater the numberof collisions in the site. According to the model by Santa Cruz et al.³²for event statistics, the combined probability for producing an eventfor protons that stop in the site and protons that traverse completelythe site is given by R/(R+2d/3), where R is the range of the proton andd is the diameter of the site. The relation between the range and theenergy can be taken from the Bragg-Kleeman rule³³, from whichR(E)=αE^(p). Taking data from NIST database³⁴ for the CSDA range ofprotons, the phenomenological function can be fitted to the ε(E) dataobtained from the MC calculations, from which can be obtained the valuesa=19.53 μm/MeV^(p) and p=1.799. From all the above, the followingfunction can be proposed to model the average energy imparted per singleevent:

$\begin{matrix}{{{\overset{¯}{ɛ}(E)} = {C_{ɛ} \cdot {{erf}( {k_{ɛ}E^{q_{ɛ}}} )} \cdot {\log( {{2{7.9}E} + b_{ɛ}} )} \cdot \frac{{1{9.5}{3 \cdot E^{({{{1.7}99} - 1})}}} + e_{ɛ}}{{1{9.5}{3 \cdot E^{{1.7}99}}} + {2{d/3}}}}},} & (9)\end{matrix}$

where C_(ε) is a proportionality constant, and b_(ε) and e_(ε) areparameters, non-physically meaningful, introduced to avoid theindetermination at E=0 and to take into account the fact thatsegment-averaged energies are being obtained, and not strictly theenergy imparted by protons of a unique energy.

As the distribution of energy imparted follows a compound Poissonprocess⁹, it is expected that its standard deviation follows a similarbehavior to its average. However, as the site diameter becomes larger,an additional effect takes place, which can be observed in FIG. 4B: aspike emerges in the low energy part of the curve for the standarddeviation. To include this effect in the model, a Gaussian term G(E) canbe subtracted from the expression in equation (9), centered at theenergy E_(c) in which the characteristic basin is observed in thosecurves. Therefore, the proposed function becomes now:

$\begin{matrix}{{{\sigma_{ɛ}(E)} = {{C_{\sigma} \cdot {{erf}( {k_{\sigma}E^{q_{\sigma}}} )} \cdot {\log( {{2{7.9}E} + b_{\sigma}} )} \cdot \frac{{19.53 \cdot E^{({{{1.7}99} - 1})}} + e_{\sigma}}{{1{9.5}{3 \cdot E^{{1.7}99}}} + {2{d/3}}}} - {C_{G}{\exp( {- {k_{G}( {E - E_{c}} )}} )}}}},} & (10)\end{matrix}$

where the Gaussian term is expressed as G(E)=C_(G) exp(−k_(G)(E−E_(c)))and incorporates three new parameters: C_(G), E_(c) and k_(G).

An example model for the average and weighted-average segment length isdescribed as follows.

For both averages of segment length, the expected behavior is asaturated function, since once the proton has a residual range largerthan the site diameter, segment length tends to the chord length and theaverages of those distributions are independent on the proton energy orrange. A good fit for this particular behavior is:

$\begin{matrix}{{\overset{\_}{s}(E)} = {\frac{2d}{3}( {{1 - {\exp( {{- k_{sF1}}{\exp( {{- k_{sF2}}E} )}} )}},} }} & (11)\end{matrix}$

where k_(sF1) and k_(sF2) are two empirical parameters and 2d/3 is thesaturation value. The same function is proposed to fit s _(D)(E) to thedata to obtain the corresponding parameters k_(sD1) and k_(sD2), butusing the weighted-average of the chord length distribution, 3d/4,instead of 2d/3 as saturation value:

$\begin{matrix}{{{\overset{\_}{s}}_{D}(E)} = {\frac{3d}{4}( {{1 - {\exp( {{- k_{sD1}}{\exp( {{- k_{sD2}}E} )}} )}},} }} & (12)\end{matrix}$

An example model for the straggling functions is described as follows.

Finally, other analytical models can be proposed for the average and thevariance—or the standard deviation—of the distribution of energyimparted per collision. As argued previously, the average energyimparted per collision should increase with the proton energy at verylow energies until reaching a maximum point in which the amount ofenergy transported out of the site by the secondary electronscompensates the increase of energy imparted due to the proton kineticenergy increase. Beyond that point, the function slowly decreases due tothe fact that the energy escape component becomes more important thanthe kinetic energy increase. A function that models well this is givenby:

$\begin{matrix}{{{\delta_{1}(E)} = {{C_{\delta_{1}}( {1 - {\exp( {k_{\delta_{1},1} \cdot E} )}} )}{\exp( {{- k_{\delta_{1},2}} \cdot E^{p_{\delta_{1}}}} )}}},} & (13)\end{matrix}$

where C_(δ) ₁ , k_(δ) ₁ _(,1), k_(δ) ₁ _(,2) and p_(δ) ₁ are the fittingparameters. Following the same reasoning than for the standard deviationof the imparted energy, σ_(δ)(E) can be modeled by fitting the samefunctional form than in equation (12) to the calculated MC data toobtain the corresponding four parameters C_(σδ), k_(σ) _(δ) _(,1), k_(σ)_(δ) _(,2), and p_(σ) _(δ) .

All the fitting to the Monte Carlo data was done with the Trust-Regionalgorithm for the method of non-linear least squares implemented inMATLAB R2018a Curve Fitting toolbox³⁵.

Example validation for other site dimensions is described as follows

Equations (9) through (12) were produced from MC simulations for sitediameters of 1 μm, 5 μm and 10 μm. In order to validate the disclosedapproach, additional MC simulations were performed for site diameters of3 μm and 7 μm and compared to the prediction provided by equations(9-12) for these site dimension.

Results of the example methods and techniques above are described asfollows. FIGS. 4A-B shows the results of the fitted functions inequations (9-10) for the average and the standard deviation of thesingle-event energy imparted for protons from 0.1 MeV to 100 MeVirradiating spherical sites with diameter of 1 μm, 5 μm and 10 μm,respectively. Residuals for these fits can be found in FIGS. 11A-B.These figures show data obtained from the MC simulations done withGeant4-DNA in liquid water for monoenergetic protons from 0.1 to 100 MeV(points) and fitted functions, equations (9-10) (lines), for sphericalsites with diameter of 1 μm, 5 μm and 10 μm, respectively. FIG. 4A showsgraphs for average single-event energy imparted to the site. FIG. 4Bshows graphs for the standard deviation of the correspondingdistributions of energy imparted. Statistical uncertainties of the MCsimulations (1σ) are shown with error bars.

As for the frequency-average segment length described by protons alongthe site, the fit of the function proposed in equation (11) to the MCcalculations is showed in FIG. 5. In this case, simulations of protonwith energies down to 10 keV can be considered, as the energy at which aproton becomes a “stopper” is below 100 keV for spherical sites of 1 μmdiameter. Residuals for this fit are shown in FIG. 12A. Results andresiduals for weighted average segment length, s _(D), equation (12),are also shown in FIG. 13 and FIG. 12B, respectively.

FIG. 5 shows graphs of calculations from MC simulations done withGeant4-DNA in liquid water for monoenergetic protons from 0.01 to 100MeV (points) and corresponding fitted function, equation (11) (lines)for the average segment length described in spherical sites withdiameters of 1 μm, 5 μm and 10 μm, respectively. Statisticaluncertainties of the MC simulation (1σ) are shown with error bars.

Once obtained the previous results, it is possible to calculate y _(D)according to the way specified in equation (2). FIG. 6 shows thecomparison between the calculated y _(D) values following this equationby using the calculations from MC simulations and the analytical fittedfunctions already shown. Differences between calculated values throughthe MC calculations and the analytical model are shown in FIG. 14.

FIG. 6 shows graphs of dose-mean lineal energy obtained from MCsimulations of monoenergetic protons in liquid water done withGeant4-DNA (points), and with the proposed analytical models (lines) byusing equation (2), for spherical sites with diameters of 1 μm, 5 μm and10 μm, respectively. Statistical uncertainties propagated from the MCsimulations (1σ) are shown with error bars.

FIG. 7 shows graphs comparing the analytical expression shown inequation (13) and the MC data resulting from Geant4-DNA simulations forthe average and the standard deviation of the distribution of impartedenergy per collision, δ₂(E). Results and residuals for the fits for eachone of these two functions are shown in FIGS. 15A-B and FIGS. 16A-Brespectively. Differences between δ₂ (E) values for MC and analyticalcalculations are shown in FIG. 17.

FIGS. 8A-B shows graphs illustrating weighted average of thedistribution of energy imparted per collision, δ₂, calculated formonoenergetic protons in liquid water from Geant4-DNA MC simulations(points), and from the analytical models described in equation (13)(lines) for spherical sites with diameters of 1 μm, 5 μm and 10 μm,respectively. Statistical uncertainties of the MC simulations (1σ) areshown with error bars.

Combining all the results shown above, for the three site diametersconsidered, it is possible to calculate the segment-averageddose-averaged restricted LET, L _(Δ,D,s), for each proton energy bymeans of equation (3). FIGS. 6A-B shows the comparison between theresults for L _(Δ,D,s) employing the data coming from the MC simulationsand from all the analytical models proposed in this work. Differencesbetween L _(Δ,D,s) values calculated from MC data and the analyticalmodels are also shown in FIGS. 8A-B.

FIG. 8A show shows graphs of segment-averaged dose-averaged restrictedLET, L _(Δ,D,s), calculated with Geant4-DNA MC simulations of protons inliquid water (points) and with the proposed analytical models (lines) byusing equation (3), for spherical site with diameters of 1 μm, 5 μm and10 μm, respectively. Statistical uncertainties propagated from the MCsimulations (1σ) are shown with error bars. FIG. 8B shows absolutedifferences between L _(Δ,D,s) calculated with the proposed analyticalmodels and with Geant4-DNA MC simulations. The point for the 10 μm siteand 0.1 MeV is not shown, being this difference equal to 22.64 keV/μm.

The difference for protons at 0.1 MeV is found to be considerably largerthan at any other energy for a site of 10 μm diameter. This is due tothe low values for s _(F) and s _(D) at that point, which implies thatsmall variations from the MC results for that quantities produce a largedeviation for L _(Δ,D,s). Leaving that point out, average differences(model calculations minus G4-DNA computations) are −0.2±0.7 keV/μm forthe 1 μm case, 0.0±0.9 keV/μm for the 5 μm case and −0.3±1.1 keV/μm forthe 10 μm case, where the uncertainties indicate the standard deviationof the differences. To look at the relevance of found differences, thesevalues can be expressed as relative to the average standard deviationfor the Monte Carlo calculations in each case. The average uncertainties(1σ) for the Monte Carlo calculations, respectively, are 0.3 keV/μm, 0.8keV/μm and 1.4 keV/μm. Consequently, all the average discrepanciesbetween the model and the MC calculations are below 1σ.

Finally, by interpolating the values of the fitting parameters obtainedin the models given by the equations (9-12), a prediction of the averageenergy imparted and L _(Δ,D,s) for intermediate site diameters can beassessed. FIGS. 9A-B shows predictions for spherical sites withdiameters of 3 μm and 7 μm, respectively, compared with the calculationsobtained from Geant4-DNA MC simulations. Residuals between the predictedresults and the MC data are shown in FIGS. 18A-B.

FIGS. 9A-B show data obtained from the Geant4-DNA MC simulations ofmonoenergetic protons in liquid water (points) and the predictionresulting from the interpolation of the fitting parameters of the models(lines) for spherical sites with diameter of 3 μm and 7 μm,respectively. FIG. 9A shows frequency-average of single-event energyimparted to the site, i; FIG. 9B shows segment-averaged dose-averagedrestricted LET L _(Δ,D,s). Statistical uncertainties propagated from theMC simulations (1σ) are shown with error bars.

The average deviation found for the L _(Δ,D,s) prediction with respectto the G4-DNA calculated points is now −1.2±1.8 keV/μm for the 3 μm siteand 1.2±3.6 keV/μm for the 7 μm, where the uncertainties refer to thestandard deviation of the differences.

Discussion of the results is provided as follows.

Three characteristics make the LET calculation presented in this workdifferent from the collision stopping power: dose-averaging, restrictionand segment-averaging. Furthermore, the calculated LET depends on twodifferent parameters: the initial proton kinetic energy and the sitedimension. Since spherical sites are used, the latter is characterizedby the site diameter. The following discussion analyzes the importanceof each characteristic of the disclosed LET calculation depending on theproton energy and the site diameter.

The dose-averaging is incorporated into the calculation by taking intoaccount the standard deviation of the different quantities involved andby using the relative variances relation leading to equation (3). It is,therefore, relevant for all proton energies and site diameters. Therestriction is performed by using a given spatial structure to score theenergy imparted and, therefore, is implicit in the distributions ofsingle-event energy imparted calculated. This spatial restriction isimportant in such situations in which secondary electrons generated bythe proton from ionizations of molecules of the medium might escape fromthe site in a non-negligible proportion. This condition is met wheneither the site is relatively small or the kinetic energy transferred tothe electrons is high enough, i.e., the proton energy is high enough tomake them escape. Therefore, one can expect to find the main differencesbetween restricted LET and stopping power at higher energies as the sitesize increases. For the same proton energy, the restricted LET value ishigher for larger sites and it tends to converge to the stopping powerwhen an infinitely large site is considered.

On the other hand, the segment-averaging process becomes important inthose cases in which the stopping power of the proton changes its valuesignificantly across the site. This happens when either the site isrelatively large or the proton energy is low enough. Secondary electronsdo not play a role in the latter case because the amount of energytransferred to them is of the order of eV, thus they are not expected totravel long distances from the proton track. It is convenient to noticethat the segment-averaging process is similar to integrate the collisionstopping power curve S(E) for protons in water along the considereddistance, i.e. the site diameter

${{\overset{\_}{L}}_{D,s} \approx {\int{{S( {E(x)} )}\frac{d{E(x)}}{dx}d{x/{\int{\frac{d{E(x)}}{dx}dx}}}}}},$

where E(x) is the energy of the proton at the point x. Therefore, thelarger is that distance, the broader and smaller results the new curve'speak with respect to collision stopping power. Additionally, because ofthe same reason, the new peak will be shifted towards higher energies.This explains why for the 1 μm diameter site such a peak is notobserved.

Additionally, the larger the site, the more the stopping power is ableto change. This means that a longer portion of the curve of stoppingpower respect to the proton energy is considered to calculate thesegment-average LET. This explains the different positions (on x-axis)and amplitude of the peaks obtained in the curves for the average andthe standard deviation of the single-event energy imparted distributionshown in FIGS. 4A-B for different site diameters. As the site becomeslarger, the peak gets higher and shifted towards higher initial protonenergies because a longer portion of the track is considered, thus aproton is able to stop and deposit all its energy in the site enteringwith higher initial energies.

The reason for the spike observed in FIG. 4B for σ_(E) at large sites(diameter of 5-10 μm) is based on the fact that distributions of energyimparted by protons at very low energies have a tailed shape at the highenergy part of the distribution, being this tail due to the eventsproduced by “stoppers”. Then, two competitive effects take place whenvery low energies are considered: as the proton energy increases, (i)the absolute energy imparted grows so that the absolute standarddeviation of the distribution does too; but (ii) the number of eventsproduced by “stoppers” decrease (their range become larger), thus thementioned tail fades away, which makes the standard deviation todecrease. For the 1 μm diameter site, the second effect occurs alreadybelow 100 keV, which is why this effect has practically no impact on thecurve. However, the larger the site, the higher the energy at which thetail disappears so that, at a certain point, the second effect starts tobe noticed, making the standard deviation decrease. Once the tailvanishes, σ_(ε) follows the dependence of f with the proton energy.

By means of equation (3), the effects of segment length variability andstraggling on the distribution of y_(s) may be disregarded in the L_(Δ,D,s) value. This is convenient because these two sources ofvariability for y_(s) carry all the dependencies with the site shapeand, as LET is a macroscopic concept, it should not depend on this siteshape. This dependency for the segment length is straightforward. In thecase of the straggling effect, as the energy imparted per collisionincludes the energy transferred by secondary electrons, the site shapeinfluences whether an electron with a given kinetic energy can leave thesite or not. According to all discussed above, therefore, L _(Δ,D,s)does not depend on the site shape but a spatial parameter given by thesite diameter d. Furthermore, the term “restricted LET” refers to aspatial restriction of the energy imparted given by the parameter A,related somehow to the site diameter d, i.e. Δ=Δ(d), and it is relevantin those cases in which secondary electrons have ranges higher orcomparable to d. On the other hand, the term “segment-averaged LET”refers to a longitudinal averaging of the energy imparted given by themean segment length, s(d). In this case, this becomes important when LETchanges significantly along that track segment.

The adequate values for the parameters Δ(d) and s(d) should be relatedto different clinical endpoints and the effect of the distribution ofenergy imparted in certain biological structures on those endpoints.Therefore, it is desirable to have a tool to calculate spatialsegment-averaged LET distributions for different sites in order to lookfor correlations between L _(Δ,D,s) values and the studied endpoints. Bystudying if these correlations exist, the quality of a restricted LETbased on a given dimension d as an indicator of some biological effectcan be assessed. In other words, LET restriction may be considered toexplain the relation between RBE and LET. This work introduces thisdegree of freedom and proposes how to handle it in order to find suchcorrelations. The concepts of unrestricted and segment-averaged LET makethis calculation expected to be different from the usual LETcalculations published in literature, which employ a concept of LETbased on the (unrestricted) electronic stopping power values for protonsin water¹²⁻¹⁵. These differences are expected to occur both at theentrance region, where the restriction effect plays a role, and at theproximities of the Bragg peak, where the segment-averaging becomesrelevant. In all cases, L _(Δ,D,s) is expected to be lower than the L_(D) calculated by those methods.

Regarding the models proposed in this work, a remarkable agreementbetween the analytical functions and the MC-based data can be found inFIGS. 3-7. This provides a double advantage: (i) calculations of L_(Δ,D,s) values can be performed much faster than MC simulations byfollowing this approach, since simple analytical functions are onlyused; and (ii) precision comparable to MC calculations, which resolvessome of the issues of the current analytical models for LETcalculations^(14,15), besides the restricted and segment-averagedcharacteristics, which, as discussed above, might link LET with clinicalendpoint. For example, performing LET calculation with this approachdoes not need to account for primary and secondary protons separately,as many of previous works do¹²⁻¹⁴. As, according to equations (6-7), onemay use the clinical beam spectral fluence ϕ(E) to convolve thepresented models and to finally obtain the corresponding LET values, itis enough with calculating the contribution from the secondaries interms of fluence instead of in terms of LET. This may represent a realadvantage when calculating lateral distribution of LET, as they onlydepend on the lateral beam fluence.

To perform clinical calculations of L _(Δ,D,s) using the proposedmethod, each domain has to have a single value of L _(Δ,D,s) which isobtained by averaging the modeled quantities respect to the spectralfluence of the beam, as shown in equations (6-7), i.e. by consideringthe distribution of entrance proton energies into the domain. Because ofthis polyenergetic averaging process, discrepancies observed in FIGS.8A-B and even in FIGS. 9A-B, which shows the prediction for sitedimensions whose data has not been employed in the model process, willprobably diffuse. This is especially likely for the high discrepanciesin low energies since in the clinical situation, relatively broadspectra are present and usually several beams converge at the samepoint. This means that the fraction of protons with low energy becomessmall, which makes the contribution of low-energy protons to theintegrals used in equations (6-7) to be irrelevant.

This work does not aim at determining the spectral fluence employed inequations (6-7), which is a separated problem to be addressed in futureworks. It should be noticed that the actual accuracy of this methodologystrongly depends on the accuracy of the determination of the spectralfluence. Essentially, this determination consists in transportingprotons in the actual media they traverse, so that either calculationsbased on MC simulations or a fluence-based kernel can be employed. As aconsequence, resolution of CT scans, for example, becomes an issue forthis spectral fluence determination. However, the intrinsic microscopicdistributions are carried by the disclosed models regardless howaccurate the spectral fluence is. On the other hand, as the disclosedmicrodosimetric models represent the way in which energy is imparted towater, the L _(Δ,D,s) calculated distributions are referred to water.This seems convenient considering that the usual way to reportdosimetric distributions in radiation therapy is dose-to-water. Thepresence of non-homogeneities, therefore, affects the proton transportproblem (at macroscopic scale) rather than the microdosimetric energydeposition.

Calculations in this work can be extended beyond the current highestproton energy simulated in Geant4-DNA, currently fixed at 100 MeV, inorder to simulate protons at the full clinical energy range. However,the contribution of protons with energy higher than 100 MeV to the LETvalue is expected to be low. Additionally, if the calculation isrestricted to the LET distribution in a tumor, those energies above 100MeV do not have much relevance. Nonetheless, an appreciable change ofbehavior in LET at 100 MeV is not expected along the rest of theclinical energy ranges, so the model presented here could beextrapolated to higher energies with probable limited error. In anycase, Monte Carlo calculations can be performed to confirm this.

Conclusions are provided as follows.

A new way of averaging LET for a particle beam, called segment-averagedLET, is added to the known operations of dose-averaging and restrictionon LET. This is done in order to incorporate the effects produced by thevariation of proton stopping power along microscopic biologicalstructures. Additionally, the quantities to calculate thissegment-averaged LET can be measured and modeled analytically frommicrodosimetric Monte Carlo simulations of monoenergetic beams. Theseanalytical models allow fast calculations of segment-averageddose-averaged restricted LET with non-significant deviations fromcalculations derived from MC simulations carried out with Geant4-DNA.From these results, a method to perform clinical calculations ofdistributions of LET over a patient based on MC calculations (e.g.,which could also be performed with a different code) can be built in astraightforward way.

The present disclosure may relate to at least the following aspects.

Aspect 1. A method comprising, consisting of, or consisting essentiallyof: determining volumetric data associated with an object of interest,wherein the volumetric data comprise a plurality of domains (e.g., thevolumetric data may be a portion of a larger set of volumetric data);determining, for each of the plurality of domains (e.g., or for aportion of the plurality of domains), a plurality of distributions forcharacterizing interactions of particles of a particle beam with thecorresponding domain, wherein the plurality of distributions comprises afirst distribution indicative of energy imparted in a correspondingdomain due to a particle traveling in the domain; determining ananalytical function based on one or more of the plurality ofdistributions; determining, based on the analytical function and foreach of the plurality of domains, first data comprising energy impartedby the particle beam to the corresponding domain; and outputting, basedon the first data, data associated with treatment of the object ofinterest by the particle beam.

Aspect 2. The method of Aspect 1, wherein the data associated withtreatment comprises one or more of (1) a three-dimensional distributionof dose or (2) a segment-averaged, restricted, dose averaged linearenergy transfer for the particle beam, and further comprising adjustinga treatment plan based on one or more of the first data or the dataassociated with the treatment.

Aspect 3. The method of any one of Aspects 1-2, wherein the first datacomprises one or more of a linear energy transfer imparted by theparticle beam to the corresponding domain or a dose imparted by theparticle beam to the corresponding domain.

Aspect 4. The method of any one of Aspects 1-3, wherein the particlebeam comprises a beam of one or more of protons, neutrons, positiveions, electrons, or alpha particles.

Aspect 5. The method of any one of Aspects 1-4, wherein determining theplurality of distributions comprises using one or more simulations togenerate the plurality of distributions.

Aspect 6. The method of any one of Aspects 1-5, wherein the plurality ofdistributions comprises a second distribution indicative of segmentlength of a particle path in the domain and a third distributionindicative of energy imparted in the domain due to a collision of aparticle.

Aspect 7. The method of Aspect 6, wherein the segment length comprises adistance a particle of the particle beam is predicted to travel afterentering the domain before coming to rest.

Aspect 8. The method of any one of Aspects 1-7, further comprisingdetermining, for a particle energy of the particle beam and based on theplurality of distributions, an energetic kernel, wherein the energetickernel comprises a first average of at least one of the plurality ofdistributions and a first variance of the at least one of the pluralityof distributions.

Aspect 9. The method of Aspect 8, further comprising performing a firstconvolution of an energy fluence of the particle beam with the firstaverage and a performing a second convolution of the energy fluence ofthe particle beam with the first variance, wherein the first data isdetermined based on a result of the first convolution and the secondconvolution.

Aspect 10. The method of any one of Aspects 1-9, wherein the dataassociated with treatment plan is based on a model that accounts for oneor more of variations of linear energy transfer in a domain, variationsof dose in a domain, variations of segment length of paths of particlesentering a domain, variations of whether particles come to rest in adomain, variations in a number of collisions of a particle in a domain,or variations in an amount of energy imparted in a collision of aparticle in a domain.

Aspect 11. The method of any one of Aspects 1-10, wherein the volumetricdata comprises one or more of geometric data associated with the objectof interest, data comprising a plurality of voxels, data associated witha cell of the object of interest, data associated with a tissue of theobject of interest, or data associated with a macroscopic structure ofthe object of interest.

Aspect 12. The method of any one of Aspects 1-11, wherein the volumetricdata comprises one or more of a geometrical model or a spatialdistribution indicative of the object of interest.

Aspect 13. The method of any one of Aspects 1-12, wherein the volumetricdata is generated based on imaging data associated with the object ofinterest.

Aspect 14. The method of any one of Aspects 1-13, wherein the domainscomprise subdivisions of a biological structure.

Aspect 15. The method of any one of Aspects 1-14, wherein one or more ofthe plurality of domains vary in one or more of shape, size, orarrangement to represent corresponding biological features of the objectof interest.

Aspect 16. The method of any one of Aspects 1-15, further comprisingdetermining data indicative of the particle beam to apply to the objectof interest.

Aspect 17. The method of Aspect 16, wherein the data indicative of theparticle beam comprises data indicative of a spatial distribution of aplurality of particles emitted from a particle emitter.

Aspect 18. The method of any one of Aspects 1-17, wherein determiningthe analytical function comprises determining the analytical functionbased on the first distribution.

Aspect 19. The method of any one of Aspects 1-18, wherein determiningthe analytical function comprises determining the analytical functionbased on fitting the one or more of the plurality of distributions tothe analytical function.

Aspect 20. The method of any one of Aspects 1-19, wherein the pluralityof domains of the volumetric data indicate a shape of one or more of acell, a nucleus of a cell, or a tissue of the object of interest.

Aspect 21. The method of any one of Aspects 1-20, wherein the pluralityof domains of the volumetric data indicate an arrangement of one or moreof a cell, a nucleus of a cell, or a tissue within the object ofinterest.

Aspect 22. A method comprising, consisting of, or consisting essentiallyof: determining data indicative of a particle beam to apply to an objectof interest; determining imaging data associated with the object ofinterest, wherein the imaging data comprises a plurality of voxels ofdata; determining, for each of a plurality of domains in a first voxelof the plurality of voxels, a plurality of distributions forcharacterizing interactions of particles of the particle beam with thecorresponding domain, wherein the plurality of distributions comprises afirst distribution indicative of energy imparted in a correspondingdomain due to a particle traveling in the domain; determining, based onthe plurality of distributions and for each of the plurality of domains,first data comprising energy imparted by the particle beam to thecorresponding domain; and outputting, based on the first data, dataassociated with treatment of the object of interest by the particlebeam.

Aspect 23. The method of Aspect 22, wherein the data associated withtreatment comprises one or more of (1) a three-dimensional distributionof dose or (2) a segment-averaged, restricted, dose averaged linearenergy transfer for the particle beam, and further comprising adjustinga treatment plan based on one or more of the first data or the dataassociated with the treatment.

Aspect 24. The method of any one of Aspects 22-23, wherein the firstdata comprises one or more of a linear energy transfer imparted by theparticle beam to the corresponding domain or a dose imparted by theparticle beam to the corresponding domain.

Aspect 25. The method of any one of Aspects 22-24, wherein the particlebeam comprises a plurality of protons at a plurality of correspondingenergies.

Aspect 26. The method of any one of Aspects 22-25, wherein determiningthe plurality of distributions comprises using one or more simulationsto generate the plurality of distributions.

Aspect 27. The method of any one of Aspects 22-26, wherein the pluralityof distributions comprise a second distribution indicative of segmentlength of a particle path in the domain and a third distributionindicative of energy imparted in the domain due to a collision of aparticle.

Aspect 28. The method of Aspect 27, wherein the segment length comprisesthe distance a particle of the particle beam is predicted to travelafter entering the domain before coming to rest.

Aspect 29. The method of any one of Aspects 22-28, further comprisingdetermining, for a proton energy of the particle beam and based on theplurality of distributions, an energetic kernel, wherein the energetickernel comprises a first average of at least one of the plurality ofdistributions and a first variance of the at least one of the pluralityof distributions.

Aspect 30. The method of Aspect 29, further comprising performing afirst convolution of an energy fluence of the particle beam with thefirst average and a performing a second convolution of the energyfluence of the particle beam with the first variance, wherein the firstdata is determined based on the results of the first convolution and thesecond convolution.

Aspect 31. The method of any one of Aspects 22-30, wherein the dataassociated with treatment plan is based on a model that accounts for oneor more of variations of linear energy transfer in a domain, variationsof dose in a domain, variations of segment length of paths of particlesentering a domain, or variations of whether particles come to rest in adomain, or variations in the number of collisions of a particle in adomain, or variations in an amount of energy imparted in a collision ofa particle in a domain.

Aspect 32. A device comprising, consisting of, or consisting essentiallyof: one or more processors; and memory storing instructions that, whenexecuted by the one or more processors, cause the device to: determinedata indicative of a particle beam to apply to an object of interest;determine imaging data associated with the object of interest, whereinthe imaging data comprises a plurality of voxels of data; determine, foreach of a plurality of domains in a first voxel of the plurality ofvoxels, a plurality of distributions for characterizing interactions ofparticles of the particle beam with the corresponding domain, whereinthe plurality of distributions comprises a first distribution indicativeof energy imparted in a corresponding domain due to a particle travelingin the domain; determine, based on the plurality of distributions andfor each of the plurality of domains, first data comprising energyimparted by the particle beam to the corresponding domain; and output,based on the first data, data associated with treatment of the object ofinterest by the particle beam.

Aspect 33. The device of Aspect 32, wherein the instructions, whenexecuted by the one or more processors, further cause the device toadjust a treatment plan based on one or more of the first data or thedata associated with the treatment.

Aspect 34. The device of any one of Aspects 32-33, wherein the firstdata comprises one or more of a linear energy transfer imparted by theparticle beam to the corresponding domain or a dose imparted by theparticle beam to the corresponding domain.

Aspect 35. The device of any one of Aspects 32-34, wherein the particlebeam comprises a plurality of protons at a plurality of correspondingenergies.

Aspect 36. The device of any one of Aspects 32-35, wherein theinstructions that, when executed by the one or more processors, causethe device to determine the plurality of distributions comprisesinstructions that, when executed by the one or more processors, causethe device to use one or more simulations to generate the plurality ofdistributions.

Aspect 37. The device of any one of Aspects 32-36, wherein the pluralityof distributions comprise a second distribution indicative of segmentlength of a particle path in the domain and a third distributionindicative of energy imparted in the domain due to a collision of aparticle.

Aspect 38. The device of Aspect 37, wherein the segment length comprisesthe distance a particle of the particle beam is predicted to travelafter entering the domain before coming to rest.

Aspect 39. The device of any one of Aspects 32-38, wherein theinstructions, when executed by the one or more processors, further causethe device to determine, for a proton energy of the particle beam andbased on the plurality of distributions, an energetic kernel, whereinthe energetic kernel comprises a first average of at least one of theplurality of distributions and a first variance of the at least one ofthe plurality of distributions.

Aspect 40. The device of Aspect 39, wherein the instructions, whenexecuted by the one or more processors, further cause the device toperform a first convolution of an energy fluence of the particle beamwith the first average and perform a second convolution of the energyfluence of the particle beam with the first variance, wherein the firstdata is determined based on the results of the first convolution and thesecond convolution.

Aspect 41. A system comprising, consisting of, or consisting essentiallyof: a particle beam generator; and at least one processorcommunicatively coupled to the particle beam generator and configuredto: determine data indicative of a particle beam to cause the particlebeam generator to apply to an object of interest; determine imaging dataassociated with the object of interest, wherein the imaging datacomprises a plurality of voxels of data; determine, for each of aplurality of domains in a first voxel of the plurality of voxels, aplurality of distributions for characterizing interactions of particlesof the particle beam with the corresponding domain, wherein theplurality of distributions comprises a first distribution indicativeenergy imparted in a corresponding domain due to a particle traveling inthe domain; determine, based on the plurality of distributions and foreach of the plurality of domains, first data comprising an energyimparted by the particle beam to the corresponding domain; and output,based on the first data, data associated with treatment of the object ofinterest by the particle beam.

Aspect 42. A device comprising, consisting of, or consisting essentiallyof: one or more processors; and a memory storing instructions that, whenexecuted by the one or more processors, cause the device to perform themethods of any one of Aspects 1-31.

Aspect 43. A non-transitory computer-readable medium storinginstructions that, when executed by one or more processors, cause adevice to perform the methods of any one of Aspects 1-31.

Aspect 44. A system comprising, consisting of, or consisting essentiallyof: a particle beam generator; and at least one processorcommunicatively coupled to the particle beam generator and configured toperform the methods of any one of Aspects 1-31.

FIG. 10 depicts a computing device that may be used in various aspects,such as in treatment planning. The computer architecture shown in FIG.10 shows a conventional server computer, workstation, desktop computer,laptop, tablet, network appliance, PDA, e-reader, digital cellularphone, or other computing node, and may be utilized to execute anyaspects of the computers described herein, such as to implement themethods described herein.

The computing device 1000 may include a baseboard, or “motherboard,”which is a printed circuit board to which a multitude of components ordevices may be connected by way of a system bus or other electricalcommunication paths. One or more central processing units (CPUs) 1004may operate in conjunction with a chipset 1006. The CPU(s) 1004 may bestandard programmable processors that perform arithmetic and logicaloperations necessary for the operation of the computing device 1000.

The CPU(s) 1004 may perform the necessary operations by transitioningfrom one discrete physical state to the next through the manipulation ofswitching elements that differentiate between and change these states.Switching elements may generally include electronic circuits thatmaintain one of two binary states, such as flip-flops, and electroniccircuits that provide an output state based on the logical combinationof the states of one or more other switching elements, such as logicgates. These basic switching elements may be combined to create morecomplex logic circuits including registers, adders-subtractors,arithmetic logic units, floating-point units, and the like.

The CPU(s) 1004 may be augmented with or replaced by other processingunits, such as GPU(s) 1005. The GPU(s) 1005 may comprise processingunits specialized for but not necessarily limited to highly parallelcomputations, such as graphics and other visualization-relatedprocessing.

A chipset 1006 may provide an interface between the CPU(s) 1004 and theremainder of the components and devices on the baseboard. The chipset1006 may provide an interface to a random access memory (RAM) 1008 usedas the main memory in the computing device 1000. The chipset 1006 mayfurther provide an interface to a computer-readable storage medium, suchas a read-only memory (ROM) 1020 or non-volatile RAM (NVRAM) (notshown), for storing basic routines that may help to start up thecomputing device 1000 and to transfer information between the variouscomponents and devices. ROM 1020 or NVRAM may also store other softwarecomponents necessary for the operation of the computing device 1000 inaccordance with the aspects described herein.

The computing device 1000 may operate in a networked environment usinglogical connections to remote computing nodes and computer systemsthrough local area network (LAN) 1016. The chipset 1006 may includefunctionality for providing network connectivity through a networkinterface controller (NIC) 1022, such as a gigabit Ethernet adapter. ANIC 1022 may be capable of connecting the computing device 1000 to othercomputing nodes over a network 1016. It should be appreciated thatmultiple NICs 1022 may be present in the computing device 1000,connecting the computing device to other types of networks and remotecomputer systems.

The computing device 1000 may be connected to a mass storage device 1028that provides non-volatile storage for the computer. The mass storagedevice 1028 may store system programs, application programs, otherprogram modules, and data, which have been described in greater detailherein. The mass storage device 1028 may be connected to the computingdevice 1000 through a storage controller 1024 connected to the chipset1006. The mass storage device 1028 may consist of one or more physicalstorage units. A storage controller 1024 may interface with the physicalstorage units through a serial attached SCSI (SAS) interface, a serialadvanced technology attachment (SATA) interface, a fiber channel (FC)interface, or other type of interface for physically connecting andtransferring data between computers and physical storage units.

The computing device 1000 may store data on a mass storage device 1028by transforming the physical state of the physical storage units toreflect the information being stored. The specific transformation of aphysical state may depend on various factors and on differentimplementations of this description. Examples of such factors mayinclude, but are not limited to, the technology used to implement thephysical storage units and whether the mass storage device 1028 ischaracterized as primary or secondary storage and the like.

For example, the computing device 1000 may store information to the massstorage device 1028 by issuing instructions through a storage controller1024 to alter the magnetic characteristics of a particular locationwithin a magnetic disk drive unit, the reflective or refractivecharacteristics of a particular location in an optical storage unit, orthe electrical characteristics of a particular capacitor, transistor, orother discrete component in a solid-state storage unit. Othertransformations of physical media are possible without departing fromthe scope and spirit of the present description, with the foregoingexamples provided only to facilitate this description. The computingdevice 1000 may further read information from the mass storage device1028 by detecting the physical states or characteristics of one or moreparticular locations within the physical storage units.

In addition to the mass storage device 1028 described above, thecomputing device 500 may have access to other computer-readable storagemedia to store and retrieve information, such as program modules, datastructures, or other data. It should be appreciated by those skilled inthe art that computer-readable storage media may be any available mediathat provides for the storage of non-transitory data and that may beaccessed by the computing device 1000.

By way of example and not limitation, computer-readable storage mediamay include volatile and non-volatile, transitory computer-readablestorage media and non-transitory computer-readable storage media, andremovable and non-removable media implemented in any method ortechnology. Computer-readable storage media includes, but is not limitedto, RAM, ROM, erasable programmable ROM (“EPROM”), electrically erasableprogrammable ROM (“EEPROM”), flash memory or other solid-state memorytechnology, compact disc ROM (“CD-ROM”), digital versatile disk (“DVD”),high definition DVD (“HD-DVD”), BLU-RAY, or other optical storage,magnetic cassettes, magnetic tape, magnetic disk storage, other magneticstorage devices, or any other medium that may be used to store thedesired information in a non-transitory fashion.

A mass storage device, such as the mass storage device 1028 depicted inFIG. 10, may store an operating system utilized to control the operationof the computing device 1000. The operating system may comprise aversion of the LINUX operating system. The operating system may comprisea version of the WINDOWS SERVER operating system from the MICROSOFTCorporation. According to further aspects, the operating system maycomprise a version of the UNIX operating system. Various mobile phoneoperating systems, such as IOS and ANDROID, may also be utilized. Itshould be appreciated that other operating systems may also be utilized.The mass storage device 1028 may store other system or applicationprograms and data utilized by the computing device 1000.

The mass storage device 1028 or other computer-readable storage mediamay also be encoded with computer-executable instructions, which, whenloaded into the computing device 1000, transforms the computing devicefrom a general-purpose computing system into a special-purpose computercapable of implementing the aspects described herein. Thesecomputer-executable instructions transform the computing device 1000 byspecifying how the CPU(s) 1004 transition between states, as describedabove. The computing device 1000 may have access to computer-readablestorage media storing computer-executable instructions, which, whenexecuted by the computing device 1000, may perform the methods describedherein.

A computing device, such as the computing device 1000 depicted in FIG.10, may also include an input/output controller 1032 for receiving andprocessing input from a number of input devices, such as a keyboard, amouse, a touchpad, a touch screen, an electronic stylus, or other typeof input device. Similarly, an input/output controller 1032 may provideoutput to a display, such as a computer monitor, a flat-panel display, adigital projector, a printer, a plotter, or other type of output device.It will be appreciated that the computing device 1000 may not includeall of the components shown in FIG. 10, may include other componentsthat are not explicitly shown in FIG. 10, or may utilize an architecturecompletely different than that shown in FIG. 10.

As described herein, a computing device may be a physical computingdevice, such as the computing device 1000 of FIG. 10. A computing nodemay also include a virtual machine host process and one or more virtualmachine instances. Computer-executable instructions may be executed bythe physical hardware of a computing device indirectly throughinterpretation and/or execution of instructions stored and executed inthe context of a virtual machine.

It is to be understood that the methods and systems are not limited tospecific methods, specific components, or to particular implementations.It is also to be understood that the terminology used herein is for thepurpose of describing particular embodiments only and is not intended tobe limiting.

As used in the specification and the appended claims, the singular forms“a,” “an,” and “the” include plural referents unless the context clearlydictates otherwise. Ranges may be expressed herein as from “about” oneparticular value, and/or to “about” another particular value. When sucha range is expressed, another embodiment includes from the oneparticular value and/or to the other particular value. Similarly, whenvalues are expressed as approximations, by use of the antecedent“about,” it will be understood that the particular value forms anotherembodiment. It will be further understood that the endpoints of each ofthe ranges are significant both in relation to the other endpoint, andindependently of the other endpoint.

“Optional” or “optionally” means that the subsequently described eventor circumstance may or may not occur, and that the description includesinstances where said event or circumstance occurs and instances where itdoes not.

Throughout the description and claims of this specification, the word“comprise” and variations of the word, such as “comprising” and“comprises,” means “including but not limited to,” and is not intendedto exclude, for example, other components, integers or steps.“Exemplary” means “an example of” and is not intended to convey anindication of a preferred or ideal embodiment. “Such as” is not used ina restrictive sense, but for explanatory purposes.

Components are described that may be used to perform the describedmethods and systems. When combinations, subsets, interactions, groups,etc., of these components are described, it is understood that whilespecific references to each of the various individual and collectivecombinations and permutations of these may not be explicitly described,each is specifically contemplated and described herein, for all methodsand systems. This applies to all aspects of this application including,but not limited to, operations in described methods. Thus, if there area variety of additional operations that may be performed it isunderstood that each of these additional operations may be performedwith any specific embodiment or combination of embodiments of thedescribed methods.

As will be appreciated by one skilled in the art, the methods andsystems may take the form of an entirely hardware embodiment, anentirely software embodiment, or an embodiment combining software andhardware aspects. Furthermore, the methods and systems may take the formof a computer program product on a computer-readable storage mediumhaving computer-readable program instructions (e.g., computer software)embodied in the storage medium. More particularly, the present methodsand systems may take the form of web-implemented computer software. Anysuitable computer-readable storage medium may be utilized including harddisks, CD-ROMs, optical storage devices, or magnetic storage devices.

Embodiments of the methods and systems are described herein withreference to block diagrams and flowchart illustrations of methods,systems, apparatuses and computer program products. It will beunderstood that each block of the block diagrams and flowchartillustrations, and combinations of blocks in the block diagrams andflowchart illustrations, respectively, may be implemented by computerprogram instructions. These computer program instructions may be loadedon a general-purpose computer, special-purpose computer, or otherprogrammable data processing apparatus to produce a machine, such thatthe instructions which execute on the computer or other programmabledata processing apparatus create a means for implementing the functionsspecified in the flowchart block or blocks.

These computer program instructions may also be stored in acomputer-readable memory that may direct a computer or otherprogrammable data processing apparatus to function in a particularmanner, such that the instructions stored in the computer-readablememory produce an article of manufacture including computer-readableinstructions for implementing the function specified in the flowchartblock or blocks. The computer program instructions may also be loadedonto a computer or other programmable data processing apparatus to causea series of operational steps to be performed on the computer or otherprogrammable apparatus to produce a computer-implemented process suchthat the instructions that execute on the computer or other programmableapparatus provide steps for implementing the functions specified in theflowchart block or blocks.

The various features and processes described above may be usedindependently of one another, or may be combined in various ways. Allpossible combinations and sub-combinations are intended to fall withinthe scope of this disclosure. In addition, certain methods or processblocks may be omitted in some implementations. The methods and processesdescribed herein are also not limited to any particular sequence, andthe blocks or states relating thereto may be performed in othersequences that are appropriate. For example, described blocks or statesmay be performed in an order other than that specifically described, ormultiple blocks or states may be combined in a single block or state.The example blocks or states may be performed in serial, in parallel, orin some other manner. Blocks or states may be added to or removed fromthe described example embodiments. The example systems and componentsdescribed herein may be configured differently than described. Forexample, elements may be added to, removed from, or rearranged comparedto the described example embodiments.

It will also be appreciated that various items are illustrated as beingstored in memory or on storage while being used, and that these items orportions thereof may be transferred between memory and other storagedevices for purposes of memory management and data integrity.Alternatively, in other embodiments, some or all of the software modulesand/or systems may execute in memory on another device and communicatewith the illustrated computing systems via inter-computer communication.Furthermore, in some embodiments, some or all of the systems and/ormodules may be implemented or provided in other ways, such as at leastpartially in firmware and/or hardware, including, but not limited to,one or more application-specific integrated circuits (“ASICs”), standardintegrated circuits, controllers (e.g., by executing appropriateinstructions, and including microcontrollers and/or embeddedcontrollers), field-programmable gate arrays (“FPGAs”), complexprogrammable logic devices (“CPLDs”), etc. Some or all of the modules,systems, and data structures may also be stored (e.g., as softwareinstructions or structured data) on a computer-readable medium, such asa hard disk, a memory, a network, or a portable media article to be readby an appropriate device or via an appropriate connection. The systems,modules, and data structures may also be transmitted as generated datasignals (e.g., as part of a carrier wave or other analog or digitalpropagated signal) on a variety of computer-readable transmission media,including wireless-based and wired/cable-based media, and may take avariety of forms (e.g., as part of a single or multiplexed analogsignal, or as multiple discrete digital packets or frames). Suchcomputer program products may also take other forms in otherembodiments. Accordingly, the present invention may be practiced withother computer system configurations.

FIGS. 11A-B show residuals for the fits of the functions in equations(9-10) to the data obtained from the MC simulations done with Geant4-DNAin liquid water for monoenergetic protons from 0.1 to 100 MeV forspherical sites with diameter of 1 μm, 5 μm and m, respectively. FIG.11A shows average single-event energy imparted to the site. FIG. 11Bshows the standard deviation of the corresponding distributions ofenergy imparted. Statistical uncertainties of the MC simulations (1σ)are shown with error bars.

FIGS. 12A-B shows residuals for the fitted functions to the calculationsfrom MC simulations done with Geant4-DNA in liquid water formonoenergetic protons from 0.01 to 100 MeV. FIG. 12A shows the averagesegment length for spherical sites with diameters of 1 μm, 5 μm and 10μm, respectively. FIG. 12B shows the weighted-average segment length forspherical sites with diameters of 1 am, 5 μm and 10 μm, respectively.Statistical uncertainties of the MC simulation (1σ) are shown with errorbars.

FIG. 13 shows calculations from MC simulations done with Geant4-DNA inliquid water for monoenergetic protons from 0.01 to 100 MeV (points) andcorresponding fitted function, equation (12) (lines) for theweighted-average segment length described in spherical sites withdiameters of 1 μm, 5 μm and 10 μm, respectively. Statisticaluncertainties of the MC simulation (1σ) are shown with error bars.

FIG. 14 shows residuals between the dose-mean lineal energy obtainedfrom MC simulations of monoenergetic protons in liquid water done withGeant4-DNA and with our proposed analytical models by using equation(2), for spherical sites with diameters of 1 μm, 5 m and 10 μm,respectively. Statistical uncertainties propagated from the MCsimulations (1σ) are shown with error bars.

FIGS. 15A-B show calculations from Geant4-DNA MC simulations formonoenergetic protons in liquid water from 0.1 MeV to 100 MeV (points)and fitted function given in equation (11) (lines). FIG. 15 A showscalculations for the average of the distribution of energy imparted percollision in spherical sites with diameters of 1 μm, 5 μm and 10 μm,respectively. FIG. 15B shows the standard deviation of the distributionof energy imparted per collision in spherical sites with diameters of 1μm, 5 μm and 10 μm, respectively. Statistical uncertainties of the MCsimulation (1σ) are shown with error bars.

FIGS. 16A-B show residuals for the fitted functions to the fromGeant4-DNA MC simulations for monoenergetic protons in liquid water from0.1 MeV to 100 MeV. FIG. 16A shows the average of the distribution ofenergy imparted per collision in spherical sites with diameters of 1 μm,5 μm and 10 μm, respectively. FIG. 16B shows the standard deviation ofthe distribution of energy imparted per collision in spherical siteswith diameters of 1 μm, 5 μm and m, respectively. Statisticaluncertainties of the MC simulation (1σ) are shown with error bars.

FIG. 17 shows residuals between the weighted average of the distributionof energy imparted per collision, δ₂, calculated for monoenergeticprotons in liquid water from Geant4-DNA MC simulations and from theanalytical models described in equation (13) for spherical sites withdiameters of 1 μm, 5 μm and 10 μm, respectively. Statisticaluncertainties of the MC simulations (1σ) are shown with error bars.

FIGS. 18 A-B show differences between data obtained from the Geant4-DNAMC simulations of monoenergetic protons in liquid water and theprediction resulting from the interpolation of the fitting parameters ofour models for spherical sites with diameter of 3 μm and 7 μm,respectively. FIG. 18A shows frequency-average of single-event energyimparted to the site, ε. FIG. 18B shows segment-averaged dose-averagedrestricted LET L _(Δ,D,s). Statistical uncertainties propagated from theMC simulations (1σ) are shown with error bars.

While the methods and systems have been described in connection withpreferred embodiments and specific examples, it is not intended that thescope be limited to the particular embodiments set forth, as theembodiments herein are intended in all respects to be illustrativerather than restrictive.

It will be apparent to those skilled in the art that variousmodifications and variations may be made without departing from thescope or spirit of the present disclosure. Other embodiments will beapparent to those skilled in the art from consideration of thespecification and practices described herein. It is intended that thespecification and example figures be considered as exemplary only, witha true scope and spirit being indicated by the following claims.

REFERENCES

-   1. Grassberger C, Paganetti H. Elevated LET components in clinical    proton beams. Phys Med Biol. 2011; 56(20):6677-6691.    doi:10.1088/0031-9155/56/20/011-   2. Paganetti H. Relative biological effectiveness (RBE) values for    proton beam therapy. Variations as a function of biological    endpoint, dose, and linear energy transfer. Phys Med Biol. 2014;    59(22):R419-R472. doi:10.1088/0031-9155/59/22/R419-   3. Belli M, Cera F, Cherubini R, et al. RBE-LET relationships for    cell inactivation and mutation induced by low energy protons in V79    cells: Further results at the LNL facility. Int J Radiat Biol. 1998;    74(4):501-509. doi:10.10100/095530098141375-   4. Peeler C R, Mirkovic D, Titt U, et al. Clinical evidence of    variable proton biological effectiveness in pediatric patients    treated for ependymoma. Radiother Oncol. 2016; 121(3):395-401.    doi:10.1016/j.radonc.2016.11.001-   5. Grün R, Friedrich T, Traneus E, Scholz M. Is the dose-averaged    LET a reliable predictor for the relative biological effectiveness?    Med Phys. 2018. doi:10.1002/mp.13347-   6. Fager M, Toma-Dasu I, Kirk M, et al. Linear Energy Transfer    Painting With Proton Therapy: A Means of Reducing Radiation Doses    With Equivalent Clinical Effectiveness. Int J Radiat Oncol. 2015;    91(5):1057-1064. doi:10.1016/j.ijrobp.2014.12.049-   7. Unkelbach J, Botas P, Giantsoudi D, Gorissen B L, Paganetti H.    Reoptimization of Intensity Modulated Proton Therapy Plans Based on    Linear Energy Transfer. Int J Radiat OncolBiolPhys. 2016;    96(5):1097-1106. doi:10.1016/j.ijrobp.2016.08.038-   8. Rossi H H, Zaider M. Microdosimetry and Its Applications.    Springer; 1996.-   9. Kellerer A M. Fundamentals of microdosimetry. In: Kase K R,    Bjarngard B E, Attix F H, eds. The Dosimetry of Ionization    Radiation. Volume I. Academic Press, Inc.; 1985:77-162.-   10. Lindborg L, Waker A. Microdosimetry. Experimental Methods and    Applications. CRC Press; 2017.-   11. Kellerer A M, Chmelevsky D. Criteria for the Applicability of    LET. Radiat Res. 1975; 63(2):226-234. doi:10.2307/3574148-   12. Wilkens J J, Oelfke U. Analytical linear energy transfer    calculations for proton therapy. Med Phys. 2003; 30(5):1006-815.    doi:10.1118/1.1567852-   13. Wilkens J J, Oelfke U. Three-dimensional LET calculations for    treatment planning of proton therapy. Z Med Phys. 2004; 14(1):41-46.    doi:10.1078/0939-3889-00191-   14. Sanchez-Parcerisa D, Cortes-Giraldo M A, Dolney D, Kondrla M,    Fager M, Carabe A. Analytical calculation of proton linear energy    transfer in voxelized geometries including secondary protons. Phys    Med Biol. 2016; 61(4):1705-1721. doi:10.1088/0031-9155/61/4/1705-   15. Marsolat F, De Marzi L, Pouzoulet F, Mazal A. Analytical linear    energy transfer model including secondary particles: Calculations    along the central axis of the proton pencil beam. Phys Med Biol.    2016; 61(2):740-757. doi:10.1088/0031-9155/61/2/740-   16. Xapsos M A. A Spatially Restricted Linear Energy Transfer    Equation. Radiat Res. 1992; 132(3):282-287.-   17. Chen J, Kellerer A M, Rossi H H. Radially restricted linear    energy transfer for high-energy protons: a new analytical approach.    Radiat Environ Biophys. 1994; 33(March):181-187.-   18. Cucinotta F A, Nikjoo H, Goodhead D T. Model for Radial    Dependence of Frequency Distributions for Energy Imparted in    Nanometer Volumes from HZE Particles Model for Radial Dependence of    Frequency Distributions for Energy Imparted in Nanometer Volumes    from HZE Particles. Radiat Res. 2000; 153(4):459-468.-   19. Rossi H H. Interpretation of Biological Response in Terms of    Microdosimetry. Ann NY Acad Sci. 1969; 161(1):260-271.    doi:10.1111/j.1749-6632.1969.tb34064.x-   20. ICRU. Report 36. Microdosimetry.; 1983.-   21. Hawkins R B. A Microdosimetric-Kinetic Model for the Effect of    Non-Poisson Distribution of Lethal Lesions on the Variation of RBE    with LET. Radiat Res. 2003; 160(1):61-69. doi:10.1667/RR3010-   22. Newhauser W D, Zhang R. The physics of proton therapy. Phys Med    Biol. 2015; 60(8):R155-R209. doi:10.1088/0031-9155/60/8/R155-   23. Kellerer A. Considerations on the random traversal of convex    bodies and solutions for general cylinders. Radiat Res. 1971;    47(2):359-376. doi:10.2307/3573243-   24. Kellerer A M. Microdosimetry and Theory of Straggling. In: Panel    on Biophysical Aspects of Radiation Quality.; 1968.-   25. Kellerer A M. Chord-Length Distributions and Related Quantities    for Spheroids. Radiat Res. 1984; 98(1):425-437.-   26. Incerti S, Baldacchino G, Bernal M, et al. The Geant4-DNA    project. Int J Model Simulation, Sci Comput. 2010; 1(2):157.    doi:10.1142/51793962310000122-   27. Incerti S, Ivanchenko A, Karamitros M, et al. Comparison of    GEANT4 very low energy cross section models with experimental data    in water. Med Phys. 2010; 37(9):4692-4708. doi:10.1118/1.3476457-   28. Bernal M A A, Bordage M C C, Brown J M C M C, et al. Track    structure modeling in liquid water: A review of the Geant4-DNA very    low energy extension of the Geant4 Monte Carlo simulation toolkit.    Phys Medica. 2015; 31(8):861-874. doi:10.1016/j.ejmp.2015.10.087-   29. Incerti S, Kyriakou I, Bernal M A, et al. Geant4-DNA example    applications for track structure simulations in liquid water: A    report from the Geant4-DNA Project. Med Phys. 2018; 45(8):e722-e739.    doi:10.1002/mp.13048-   30. Bethe H. Zur Theorie des Durchgangs schneller    Korpuskularstrahlen durch Materie. Ann Phys. 1930; 397(3):325-400.    doi:10.1002/andp.19303970303-   31. ICRU. Report 73. Stopping of ions heavier than Helium. J Int    Comm Radiat Units Meas. 2005; 5(1):iii-viii.    doi:10.1093/jicru/ndi002-   32. Santa Cruz G A, Palmer M R, Matatagui E, Zamenhof R G. A    theoretical model for event statistics in microdosimetry. I: Uniform    distribution of heavy ion tracks. Med Phys. 2001; 28(6):988-996.    doi:10.1118/1.1376439-   33. Ulmer W, Matsinos E. Theoretical methods for the calculation of    Bragg curves and 3D distributions of proton beams. Eur Phys J Spec    Top. 2010; 190(1):1-81. doi:10.1140/epjst/e2010-01335-7-   34. Berger M J, Coursey J S, Zucker M A, Chang J. Stopping-Power;    Range Tables for Electrons, Protons, and Helium Ions | NIST. NIST    Standard Reference Database 124.    https://www.nist.gov/pml/stopping-power-range-tables-electrons-protons-and-helium-ions.    Published 2017. Accessed Dec. 19, 2018.-   35. The MathWorks Inc. Matlab R2018a and Curve Fitting toolbox.    2018.

Additional Examples

An analytical microdosimetric model for radioimmunotherapeutic alphaemitters is described below.

The disclosed provides a methodology to analytically determinemicrodosimetric quantities in radioimmunotherapy and targetedradiotherapy with alpha particles.

Methods and material: Monte Carlo simulations (MC) using Geant4-DNAtoolkit, that provides physics at the microscopic level, are performedfor monoenergetic alpha particles traversing spherical sites withdiameters of 1 μm, 5 μm and 10 μm. An analytical function is fittedagainst the data in each case to model the energy imparted by eachmonoenergetic particle into the site, as well as the variance of thedistribution of energy imparted. Those models allow to obtain the meanand dose-mean values of specific energy (z) and lineal energy (y) forpolyenergetic arrangements of alpha particles. The energetic spectrum isestimated by considering the distance that each particle needs to travelto reach the target. We apply this methodology to a simple case inradioimmunotherapy: a spherical cell that has its membrane uniformlycovered by ²¹¹At, an alpha emitter, with a spherical target representingthe nucleus, placed at the center of the cell. We compare the results ofour analytical method with calculations with Geant4-DNA of this specificsetup for three nucleus sizes corresponding to our three functions.

Results: For nuclei of 1 μm and 5 μm of diameters, all mean and dosemean quantities for y and z were in an agreement within 4% to Geant4-DNAcalculations. This agreement improves to about 1% for dose-mean linealenergy and dose-mean specific energy. For the 10-μm diameter case,discrepancies scale to about 9% for mean values and 3% for dose-meanvalues. Dose-mean values are within Geant4-DNA uncertainties in allcases.

Conclusions: Our method provides accurate analytical calculations ofdose-mean quantities that may be further employed to characterizeradiobiological effectiveness of targeted radiotherapy. The spatialdistribution of sources and targets is required to calculatemicrodosimetric relevant quantities.

Radioimmunotherapy is a modality for cancer treatment based on the useof radiolabeled antibodies highly affine to antigens particularlyexpressed in tumor cell environments (R1). Specifically, alpha emissionsenhance the radiobiological effectiveness of a given dose due to theirhigher linear energy transfer (LET) in contrast with beta emissions(R2). Among others, some radioimmunotherapeutic alpha particle-emittingradionuclides are ²¹¹At, ²¹²Bi, ²²⁵Ac, ²²³Ra or ²⁵⁵Fm (R3), all of themwith emissions of energies ranging from 5.9 MeV to 9.0 MeV and diverseconstraints such as availability, radiation safety or conjugatestability (R4). The range in liquid water for alpha particlescorresponding to these energies spans approximately from 50 μm to 100μm, whereas the range for typical radioimmunotherapeutic beta emissionsis of the order of millimeters (R5). Thus, applications ofalpha-emitters in radioimmunotherapy generally exploit this highlylocalized energy deposition, for example, to treat bone metastases with²²³Ra dichloride (R6) or to target cancer cells in transit in thevascular and lymphatic systems (R7).

Generally, the main quantity in radiotherapy to predict the biologicaleffect is the dose absorbed by a tissue or other macroscopic structure.Thus, as long as this cell-wise distribution is relatively homogeneous,the mean of a distribution of doses absorbed by the cells compoundingthat tissue can be taken as a predictor for the effect. However, such ashort range for radioimmunotherapeutic alpha-emitters produces a highlyinhomogeneous cell-wise distribution of dose, so that the macroscopicdose is no longer valid to characterize the biological effect (R8).Therefore, a description of dose at a microscopic level is required.Microdosimetry provides a framework that fits this context (R9, R10).

A series of studies has recently been published (R11-R13) with thepurpose of modeling the microdosimetric patterns of energy depositionfor external proton radiotherapy. In this work, we apply the sameprinciples to alpha particles in order to produce analytical functionsthat characterize the microdosimetric distributions for monoenergeticbeams of alpha particles. Although other analytical approaches have beenproposed to implement microdosimetry for the use of alpha-emitters inradioimmunotherapy (R14-R16), the methods developed throughout thisseries are here adapted to provide a simpler and faster method tocalculate microdosimetric quantities in this application. Thisdetermination enables biophysical models to be used inradioimmunotherapy (R17-R19).

A particular set of conditions of microdosimetry for radioimmunotherapyis that the spatial distributions of (i) the emitting alpha sources, and(ii) the target cells in a tissue or other macroscopic structure need tobe known due to the short range of the emitted alpha particles. On theone hand, depending on the specific application, the spatialdistribution of emitters can be imaged by autoradiography (R20),histologically (R21) or assumed as homogeneous on the cell membrane orthe cellular media (R22-R25). On the other hand, theapplication-dependent distribution of targets throughout tissues, andeven inside of a cell (R26)—can also be imaged (R27) or modeled (R28).In this work, we present our analytical model for alpha particles andrevisit our general procedure to carry out microdosimetric calculations,applied to immunotherapy using simple geometries for aspatially-homogeneous distribution of sources.

Methods and Materials

Principles of Microdosimetry

At the microscopic scale, the nature of the interaction betweenradiation and matter is stochastic. This means that, for a given setup,the amount of energy imparted E to a certain micrometer-sized volume,i.e. a site, varies for each particle track, i.e. event. The quotientbetween ε and the mass of the site defines the specific energy, z≡ε/m.The specific energy for a single event is represented by z₁≡ε_(s)/m,where εs is the energy imparted to the site in a single event (R29).After a number of events a distribution for z₁, f₁(z), can be obtained.Also, if the same experiment is repeated several times with a fixednumber of events, i.e. of particles, another distribution for z, f(z),can be obtained. The mean of this distribution coincides with themacroscopic concept of dose: D=z.

Lineal energy (y) is a different microdosimetric quantity defined as thequotient between energy imparted to the site in a single event and themean chord length l of the particles track within the site (R30):y≡ε_(s)/l. This traditional definition is valid as long as the range ofthe considered particles is very long with respect to the site typicaldimensions (R31). However, in radioimmunotherapy, internalalpha-emitters produce particles with ranges of the order of thebiologically relevant sites. To overcome this limitation, the concept ofsegment length s and segment-based lineal energy y_(s)≡ε/s have beenintroduced (R12) and single-event distributions for y_(s), f(y_(s)), canbe obtained. The descriptions based on specific energy and lineal energyare related by:

$\begin{matrix}{z_{1} = {\frac{\overset{\_}{s}}{m}y_{s}}} & ({A1})\end{matrix}$

For the sake of simplicity, we call y≡y_(s) hereinafter. According tothe Theory of Dual Radiation Action (TDRA), sub-lethal lesions producedby independent tracks can interact to form new lethal lesions (R32). Inother words: not only does the biological effect depend on the firstmomentum of the distribution of ε_(s) (equivalently, z₁ or y_(s)) butalso it depends on the second momentum. Therefore, the weighted averages(also called dose-weighted means) of these distributions, definedgenerically as x _(D)=x² /x where x may represent ε_(s), z₁ or y_(s),become relevant to characterize the effect produced by a givenradiation. Thus, the calculation of the quantity z _(1,D) in typicalradioimmunotherapy conditions or, equivalently according to equation(A1), y_(D), is the aim of this work.

Framework for Microdosimetric Calculations

The dose-weighted mean single-event specific energy, z _(1,D), can beexpressed, according to the definition of z₁, as z _(1,D)=ε _(s,D)/m,where ε _(s,D) is the dose-weighted mean energy imparted per event tothe site. In turn, ε _(s,D) can be decomposed in terms of the meanenergy imparted per event ε _(s) and its variance, σ_(ε) _(s) ², bymeans of the identity (R11)

$\begin{matrix}{{\overset{¯}{ɛ}}_{s,D} = {{\overset{¯}{ɛ}}_{s}( {1 + \frac{\sigma_{ɛ_{s}}^{2}}{{\overset{¯}{ɛ}}_{s}^{2}}} )}} & ({A2})\end{matrix}$

Let us assume that the energy density distribution of the beam cominginto the site, ϕ_(E)(E)≡dϕ(E)/dE, is known, where E is the particlekinetic energy and ϕ(E) is the integral fluence of the beam, i.e. thenumber of particles with energy lower than or equal to E. Then, if themean energy imparted per event and the variance of f(ε_(s)) for amonoenergetic beam with energy E is given by ε _(s)(E) and σ_(ε) _(s)²(E), respectively, the corresponding quantities for a polyenergeticbeam is given by (R12)

$\begin{matrix}{{\overset{\_}{ɛ}}_{s} = {\frac{\int{{\phi_{E}(E)}{{\overset{\_}{ɛ}}_{s}(E)}dE}}{\int{{\phi_{E}(E)}dE}}\mspace{14mu}{and}}} & ({A3}) \\{\sigma_{ɛ_{s}}^{2} = {\frac{\int{{\phi_{E}(E)}{\sigma_{ɛ_{s}}^{2}(E)}dE}}{\int{{\phi_{E}(E)}dE}} + \frac{\int{{\phi_{E}(E)}( {{{\overset{\_}{ɛ}}_{s}(E)} - {\overset{\_}{ɛ}}_{s}} )^{2}dE}}{\int{{\phi_{E}(E)}dE}}}} & ({A4})\end{matrix}$

respectively. Consequently, if the energy distribution, or spectralfluence, of the particles arriving at the site is provided, z _(1,D) canbe calculated combining equations (A2)-(A4) through the analyticalenergy-dependent functions ε _(s)(E) and σ_(ε) _(s) ²(E).

Analytical Model for Monoenergetic Functions

In a similar way as shown in Bertolet et al. (R12),semi-phenomenological analytical functions for ε _(s)(E) and σ_(ε) _(s)²(E) are proposed based on the following principles: (1) both the meanenergy imparted per event and its standard deviation are proportional tothe stopping power of the particle in water and inversely proportionalto the energy; (2) as energy increases, secondary electrons increasetheir range and are more likely to escape the site, reducing the energydeposited by the track within the site; and (3) for very low energies orshort ranges, the probability for a particle going towards the peripheryof the site to reach it decreases due to the sphere curvature. Afunction that fulfills these principles is given by (R12)

$\begin{matrix}{{{\overset{\_}{ɛ}}_{s}(E)} = {er{{f( {k_{s}E^{q_{s}}} )} \cdot {\log( {{a \cdot E} + b_{s}} )} \cdot \frac{{\alpha \cdot E^{p - 1}} + e_{s}}{{\alpha \cdot E^{p}} + {4{r_{s}/3}}}}}} & ({A5})\end{matrix}$

where C_(s), k_(s), q_(s), b_(s) and e_(s) are fitting parametersdependent on the site dimension, r_(s) is the radius of the site and a,α and p are physical parameters particle-dependent. These parameters foralpha particles take the values a=7.017 MeV⁻¹, α=1.946 μm/MeV^(p) andp=1.752 following the same reasoning as in Bertolet et al. (R12) foralpha particles instead of protons. The standard deviation of f(ε_(s))follows a similar behavior (R12) so that we use the function

$\begin{matrix}{{\sigma_{ɛ}(E)} = {C_{s}^{\prime} \cdot {{erf}( {k_{s}^{\prime}E^{q\prime_{s}}} )} \cdot {\log( {{a \cdot E} + b_{s}^{\prime}} )} \cdot \frac{{\alpha \cdot E^{p - {f\prime}_{s}}} + e_{s}^{\prime}}{{\alpha \cdot E^{p}} + {4{r_{s}/3}}}}} & ({A6})\end{matrix}$

where a, α, p and d are the same as in equation (A5) while C′_(s),k′_(s), q′_(s), b′_(s), f′_(s) and e′_(s) are fitting parametersdifferent from equation (A5). Here we introduce a generalization of thearguments and equations shown in our previous work for protons to modelthe microdosimetry of alpha particles.

Monte Carlo Simulations for Monoenergetic Alpha-Particle Beams

Monoenergetic alpha-particle beams were simulated in liquid water withthe Monte Carlo (MC) code Geant4-DNA (R33-R36) by employing anapplication similar to that shown in previous works (R11,R12): alphaparticles with energies up to 50 MeV come from a point source to travelthrough a liquid water box in which the position of a spherical site issampled with uniform probability within the box. To ensure the conditionof charged particle equilibrium, margins are added upstream anddownstream with thickness equal to the maximum range of the secondaryelectrons generated by the electronic collisions of the alpha particleswith water. More details about the geometry and setup of thesesimulations can be found in Bertolet et al. (R12). These simulationswere repeated for spherical sites of 1 μm, 5 μm and 10 μm diameter.

Particle-Alpha Spectrum for Radioimmunotherapy

Given a spatial arrangement of radiolabeled antibodies isotropicallyemitting alpha particles and a spherical target placed at a certainposition (points on a cellular membrane in this work), a distribution ofdistances d between sources and the site can be calculated. Note that,to reproduce the conditions in which analytical functions from equations(A5) and (A6) were obtained, these distances need to be calculated withrespect to the plane tangential to the sphere and perpendicular to theparticle track. Here we consider the range of the particle as that givenby continuous slowing down approximation (CSDA). As radionuclides emitalpha particles with fixed energies E_(i), the spectrum ofalpha-particles arriving at the site can be calculated as

$\begin{matrix}{{\phi_{E}(E)} = {\phi_{E}( {E( {{R( E_{i} )} - d} )} )}} & ({A7})\end{matrix}$

where R(E_(i)) is the range for alpha particles with energy E_(i), d isa variable that represents the distribution of distances between sourcesand site, and E(x) is the energy of an alpha particle corresponding to arange x. Given a certain point along the particle trajectory, we callresidual range to the CSDA range that the particle still has at thatpoint. In other words, the residual range at a point is the initialrange minus the already travelled distance d: R(E_(i))−d. Finally, theASTAR database from the NIST (R37) can be used to obtain the energycorresponding to that CSDA residual range, (E(R(E_(i))−d)).

Application to a Simple Geometry

If the spatial distributions of sources and the target position areknown, i.e. f(d) is known, z _(1,D) can be calculated by using equations(A1) to (A7). Here, we illustrate our method and check our results for asimple geometrical configuration previously used in other works (R14): aspherical target inside a spherical cell with an homogeneousdistribution of sources throughout the cellular membrane.

FIGS. 19A-B show the configuration for the first case. A spherical cellwith radius r_(c) contains the target, another sphere with radius r_(n)at its center, and radioactive sources are uniformly distributed acrossthe cellular membrane surface.

FIG. 19A shows radionuclides uniformly distributed throughout themembrane of a spherical cell with radius rc. The nucleus (inner spherewith radius rn) is placed at the center of the cell and it is consideredhere as the target. As radionuclides irradiate isotropically, from eachpoint on the membrane a cone of lengths arrives at the target. FIG. 19Bshows the distance d for each particle coming out of a membrane point isthe distance to a plane perpendicular to the particle track andtangential to the nucleus. As the schematic shows, for an angle θbetween the particle track and the cell radius, the distance to thisplane is given by d=r_(c) cos θ−r_(n).

The distance d that an alpha particle travels until reaching the planetangential to the target and normal to its track is given by d=r_(c) cosθ−r_(n), where θ is the angle with respect to surface normal of themembrane (see FIG. 19B). As each source emits isotropically, the numberof alpha particles for each angle is proportional to the area dA=2πr drof a ring with radius r and thickness dr perpendicular to a cell radius,so that r=rc tan θ. Then, the angular distribution of tracks from eachmembrane point is given by

$\begin{matrix}{{{g(\theta)}d\theta} \propto {2\pi r_{c}^{2}\frac{\sin(\theta)}{\cos^{3}(\theta)}d\theta}} & ({A8})\end{matrix}$

which can be expressed in terms of the distance d using the previousrelation as

$\begin{matrix}{{{f(d)}dd} = {{N \cdot \frac{2\pi}{( {d + r_{n}} )^{3}}}dd}} & ({A9})\end{matrix}$

where N is a normalization factor so that ∫_(d) _(min) ^(d) ^(max)f(d)dd=1/N, d_(max) is the maximum distance, which corresponds to thecase θ=0, in which d_(max)=r_(c)−r_(n), and d_(min) is the minimumdistance, that is found when cos θ=√{square root over (1−r_(n) ²/r_(c)²)}, so that d_(min)=√{square root over (r_(c) ²−−r_(n) ²)}−r_(n). Thisintegration yields N=(r_(c) ²−r_(n) ²)/πr_(n) ². On the other hand, asall the points on the membrane are in a symmetrical position, thedistribution of distances for all the sources is given by

$\begin{matrix}{{{f(d)}dd} = {\frac{2( {r_{c}^{2} - r_{n}^{2}} )}{{r_{n}^{2}( {d + r_{n}} )}^{3}}dd}} & ({A10})\end{matrix}$

If we consider, for example, ²¹¹At as the emitting radionuclides, thespectrum of alpha particles is composed of 42% of 5.87-MeV particles and58% of 7.45-MeV particles (R38), or, respectively 48.0 μm and 69.9 μm inrange. Then, we can obtain the distribution of residual ranges bysubtracting the distribution of distances d to these ranges and,finally, the energy density distribution ϕ_(E)(E) by converting theseresidual ranges into energy as pointed out in equation (A7). With afixed radius of 7.5 μm for the cell, we calculate the analytical resultsfor three different nucleus radii: 0.5 μm, 2.5 μm and 5.0 μm. Forexample, for the 2.5 μm case, the maximum distance becomes 5 μm and theminimum 4.57 μm. The resulting spectrum for this case is shown in FIG.20.

FIG. 20 shows a spectrum ϕ_(E)(E) of alpha particles coming to thespherical target of radius r_(n)=2.5 μm shown in FIG. 19B for a uniformactivity of ²¹¹At distributed around the membrane of a spherical cell ofradius r_(c)=7.5 μm.

Monte Carlo Simulations of the Radioimmunotherapeutic Example

In order to benchmark our analytical calculation, we have reproduced theproblem described in FIGS. 19A-B with the MC toolkit Geant4 (R39-R41)(v10.5) and the track-structure physics built in the package Geant4-DNA.In our Geant4 application, a spherical surface with 211At sources isuniformly distributed with a spherical target volume at its center withthe same radii as specified above. For each radius size we simulate225,000 alpha tracks, considering a uniform angular distribution butonly directed to the solid angle subtended by the nucleus (as shown inFIG. 19A) plus a margin given by the maximum range of secondaryelectrons generated by alpha particles with the considered energies.This margin is added to include energy depositions by secondaryelectrons entering the nucleus when the primary track passes near itssurface, and, for example, for the case with radius of 2.5 μm, it rangesbetween 0.34 μm and 0.51 μm. The cutoff production for secondaryparticles were set to 0.5 μm range. Geant4-DNA models were activatedbelow 1 MeV as kinetic energy of electrons and 10 MeV for protons andalpha particles. Particularly we used for electrons Champion model forelastic scattering and Born model for ionization and excitationprocesses. As for protons and alphas, we used: (a) for kinetic energiesbelow 0.5 MeV: Bragg ion gas model for elastic processes, Born model forionization and Miller-Green model for excitation; and (b) for energiesabove 0.5 MeV: Bethe-Block ion gas model for elastic processes and Bornmodel for ionization and excitation. This selection of Geant4-DNA modelsis based on the Geant4 official example extended/medical dnamicrodosimetry.

As this geometrical configuration does not provide a uniform randomnessfor the intersection of particle tracks and site (R42), we obtain y forthe simulations by collecting the distribution of segment length and itsmean value so that y_(s)≡y=ε_(s)/s. Then, mean values for z₁ arecorrected in our analytical calculation according to equation (A1) byusing the obtained s. For the MC simulations, z₁ is simply obtained asthe scored quotient ε_(s)/m, being m the mass within the site.

Results

Models for Alpha Particles

FIGS. 21A-B shows a comparison among the data gathered frommonoenergetic simulations of alpha particles in Geant4-DNA and theproposed analytical functions in equations (A5) and (A6) for sites withdiameters of 1 μm, 5 μm and 10 μm, respectively.

FIGS. 21A-B show example analytical functions (solid line) fitted to thedata (dots) obtained from monoenergetic Geant4-DNA simulations. FIG. 21Ashows example analytical functions for mean energy imparted to sphericalsites with diameters of 1 μm (top), 5 m (center) and 10 μm (bottom),respectively. FIG. 21B shows example analytic functions for variance ofthe energy imparted to spherical sites with diameters of 1 μm (top), 5μm (center) and 10 μm (bottom), respectively. Statistical uncertaintiesfrom simulations (1σ) are shown with error bars.

Results for the Analytical Calculation and MC Simulation

We calculate ε _(s) and σ_(ε) ² by using the calculated spectrumϕ_(E)(E), as shown for instance in FIG. 20 for the nucleus of 5 μmin-diameter, and the modelled functions shown in FIGS. 21A-B. Then, weobtain the mean values for the single-event distribution of specificenergy z₁ and lineal energy y as explained above. For the threedifferent nucleus radii, the mean segment length scored in theGeant4-DNA simulations were (a) 0.660±0.001 μm for r_(n)=0.5 μm; (b)3.155±0.003 μm for r_(n)=2.5 μm; and (c) 5.50±0.09 μm for r_(n)=5.0 μm.The comparison between analytical results and Geant4-DNA simulations areshown in FIGS. 22A-B. Statistical uncertainties for the Geant4-DNA wereestimated according to the method explained in the Appendix of (R11).Note that while y_(F) only depends on the mean energy imparted to thesite, y_(D) also depends on the variance of energy imparted to the site.As the site becomes larger by a factor x, both mean energy and segmentlength increase roughly by the same factor x, but variance of energyimparted increases as x². This implies that the variability on y_(D)becomes larger the larger the site is.

A summary of the observed differences for the mean values of y and z₁calculated analytically and through MC computations is provided inTable 1. For nucleus radii corresponding to 0.5 μm and 2.5 μm, alldiscrepancies for the considered quantities are below 4%, withespecially good agreements for dose-mean quantities. The discrepanciesare larger for the 5.0 μm-case, which indicates a better performance ofthe functions for 1 μm and 5 μm site diameters compared with those for10 μm site diameter, shown in FIG. 20.

TABLE 1 y _(F) (keV/μm) y _(D) (keV/μm) z _(1, F) (Gy) z _(1, D) (Gy)r_(n) = 0.5 μm  1.3 ± 0.2 [2.2%]  −1 ± 4 [−1.2%]  0.27 ± 0.03 [2.2%]−0.2 ± 0.9 [−1.1%] r_(n) = 2.5 μm 2.22 ± 0.11 [3.3%]  0 ± 11 [−0.9%]0.017 ± 0.001 [3.2%] 0.00 ± 0.09 [−0.9%]  r_(n) = 5.0 μm 5.69 ± 0.11[8.5%] −3 ± 15 [−3.2%] 0.010 ± 0.001 [8.9%]  0.0 ± 0.4 [−3.3%]

Discussion

We have adapted our analytical functions derived for the microdosimetricbehavior of protons in water in external radiation therapy treatments(R12) to the alpha particle cases by modifying some of the parametersrelated to the physical characteristics of each particle. The samefunctions as used for protons can be fitted to the microdosimetric MCresults for alpha particles as the underlying physical processes areessentially the same.

FIGS. 22A-B show results from analytical calculations and Geant4-DNAsimulations for a spherical cell of radius 7.5 μm and variable nucleusradius. FIG. 22A shows analytical calculations and Geant4-DNAsimulations for a spherical cell of radius 7.5 μm and variable nucleusradius y _(F) and y _(D). FIG. 22B shows analytical calculations andGeant4-DNA simulations for a spherical cell of radius 7.5 μm andvariable nucleus radius z_(1,F) and z_(1,D). Statistical uncertaintiesfor Geant4-DNA simulations show 1σ and are obtained as explained in theAppendix of Bertolet et al. 2019 (R11). Note that the scale for z1values is logarithmic in order to facilitate the visualization of theresults.

Although, as said above, the descriptions of microdosimetry based on yor z₁ are essentially equivalent, we have presented the results for bothof them. As the formulation of the TDRA is based on the distribution ofspecific energy, FIG. 22A shows a practical property of the linealenergy-based description: y _(D) is relatively insensitive to changes ofsite geometries. This refers to changes of site dimensions keeping agiven site shape, as shown by the results from Geant4-DNA simulations,but also to changes in the randomness of the segment lengthdistribution. On the one hand, as y _(D) seems relatively constantregardless the site size according to FIG. 22A, it is potentiallypossible to just use a single model from the three presented in FIGS.21A-B to calculate y _(D) corresponding to any other site size.Nonetheless, further study is needed due to the large uncertainties fory _(D) values. On the other hand, our analytical models are based onsituations of uniform randomness, i.e., alpha particles traverse thesite with a distribution of lines placed uniformly in space. However,the setup shown in FIGS. 19A-B does not contain the same assumption asall lines traversing the site come from the same point, which isreflected on the obtained values for the mean segment length, differingfrom the s=4r_(d)/3 expected in uniform randomness situations. Theagreement between our analytical model and Geant4-DNA points out that y_(D), as it is specific to the unit length, are insensitive even tochanges of segment length distributions, at least in the casesconsidered here. In other words, y _(D) can be calculated from ourspherical models regardless the actual segment length randomness. Then,z _(D) can be calculated for the specific setup considered by usingequation (A1) with the actual s obtained. Consequently, two differentdistributions need to be estimated to be able to correctly apply ourmethod: (i) the distribution of distances between sources and targets inorder to estimate the energetic spectrum; and (ii) the distribution ofsegment lengths in the target to correctly transform thegeometrical-insensitive means of y into radiobiologically relevant meansof z.

Another factor potentially affecting the results of our analyticalmethod is the dependence between the incoming spectrum and thedistribution of segment length. According to FIG. 19B, distance is afunction of the angle θ subtended by the cell radius and the directionof the alpha particle. Distance to the plane tangential to the nucleusfor a particle going toward the center (θ=0) is longer than distancesfor particles toward the periphery. Therefore, at the inner regionsalpha particles are less energetic than at the outer region whenarriving at such tangential plane: i.e. particles with longer segmentlength have also less energy whereas in the analytical model it isassumed a uniform distribution of segment length for each spectralcomponent. However, this does not seem to have a major impact ondose-mean values, probably because the width of the differences in theenergy imparted by different components of the incoming spectrum issmall.

Further works based on this analytical calculation might incorporateradiobiological models for the direct damage to the cell to convertphysical dose or microdosimetric quantities into actual biologicalpredictions, which remains as one of the challenges in targetedradiotherapy (R43). In this sense, any model would require thedetermination of radiobiological parameters, such as α and β from theLinear-Quadratic (LQ) model. The determination of these parametersusually carries large experimental uncertainty and considerablevariability among experiments can be observed in literature (R44).Therefore, the discrepancies obtained in this work for z1,D are probablysmall compared with that source of uncertainty.

CONCLUSIONS

An adaptation for alpha particles of the previously presentedmethodology to calculate microdosimetric quantities of biophysicalinterest allows a simple and direct application of microdosimetry to thefield of radioimmunotherapy. In order to apply our analytical approach,it is necessary to determine the spatial distribution of sources andtheir distances to the considered targets in order to compute theenergetic spectrum and the segment length distribution. We show a goodagreement between the values obtained from analytical microdosimetriccalculations and those obtained from Geant4-DNA for the case of uniformactivity distributed upon a spherical surface with a spherical targetinside. This methodology may be further applied to quantification ofbiological effectiveness of direct damage in radioimmunotherapy byemploying models with basis on microdosimetry and more biologicallyrelevant geometries.

ADDITIONAL REFERENCES

-   R1. Knox S J, Levy R, Miller R A, Uhland W, Schiele J, Ruehl W, et    al. Determinants of the Antitumor Effect of Radiolabeled Monoclonal    Antibodies. Cancer Res. 1990; 50(16):4935-40.-   R2. Langmuir V K, Fowler J F, Knox S J, Wessels B W, Sutherland R M,    Wong J Y C. Radiobiology of radiolabeled antibody therapy as applied    to tumor dosimetry. Med Phys. 1992; 20(2):601-9.-   R3. McDevitt M R, Sgouros G, Finn R D, Humm J L, Jurcic J G, Larson    S M, et al. Radioimmunotherapy with alpha-emitting nuclides. Eur J    Nucl Med. 1998; 25(9):1341-51.-   R4. Imam S K. Advancements in cancer therapy with alpha-emitters: A    review. Int J Radiat Oncol Biol Phys. 2001; 51(1):271-8.-   R5. Zalutsky M R, Bigner D D. Radioimmunotherapy with α-particle    emitting radioimmunoconjugates. Acta Oncol (Madr). 1996;    35(3):373-9.-   R6. Jadvar H. Targeted radionuclide therapy: An evolution toward    precision cancer treatment. Am J Roentgenol. 2017; 209(2):277-88.-   R7. Allen B, Huang C-Y, Clarke R. Targeted alpha anticancer    therapies: update and future prospects. Biol Targets Ther. 2014;    255.-   R8. Hofmann W, Li W B, Friedland W, Miller B W, Madas B, Bardiès M,    et al. Internal microdosimetry of alpha-emitting radionuclides    [Internet]. Radiation and Environmental Biophysics. Springer Berlin    Heidelberg; 2019. Available from:    https://doi.org/10.1007/s00411-019-00826-w-   R9. Kellerer A M. Fundamentals of microdosimetry. In: Kase K R,    Bjarngard B E, Attix F H, editors. The Dosimetry of Ionization    Radiation Volume I. Academic Press, Inc.; 1985. p. 77-162.-   R10. Rossi H H, Zaider M. Microdosimetry and Its Applications.    Springer; 1996.-   R11. Bertolet A, Baratto-Roldin A, Barbieri S, Baiocco G, Carabe A,    Cortés-Giraldo M A. Dose-averaged LET calculation for proton track    segments using microdosimetric Monte Carlo simulations. Med Phys    [Internet]. 2019 Jul. 12; 46(9):4184-92. Available from:    https://onlinelibrary.wiley.com/doi/abs/10.1002/mp.13643-   R12. Bertolet A, Baratto-Roldin A, Cortes-Giraldo M A,    Carabe-Fernandez A. Segment-averaged LET concept and analytical    calculation from microdosimetric quantities in proton radiation    therapy. Med Phys [Internet]. 2019; 46(9):4204-14. Available from:    http://doi.wiley.com/10.1002/mp.13673-   R13. Bertolet A, Cortes-Giraldo M A, Souris K, Cohilis M,    Carabe-Fernandez A. Calculation of clinical dose distributions in    proton therapy from microdosimetry. Med Phys. 2019 Oct. 11;    46(12):5816-23.-   R14. Stinchcomb T G, Roeske J C. Analytic microdosimetry for    radioimmunotherapeutic alpha emitters. Med Phys. 1992;    19(6):1385-93.-   R15. Stinchcomb T G, Soyland C, Hassfjell S P, Westman J, Wang S J,    Whitlock J L, et al. Binary methods for the microdosimetric analysis    of cell survival data from alpha-particle irradiation. Cancer    Biother Radiopharm. 2003; 18(3):481-7.-   R16. Chouin N, Bitar A, Lisbona A, Chérel M, Davodeau F, Barbet J,    et al. Implementation of a microdosimetric model for    radioimmunotherapeutic alpha emitters. Cancer Biother Radiopharm.    2007; 22(3):387-92.-   R17. Roeske J C, Stinchcomb T G. Tumor Control Probability Model for    Alpha-Particle-Emitting Radionuclides. Radiat Res. 2000;    153(1):16-22.-   R18. Roeske J C, Stinchcomb T G. The average number of    alpha-particle hits to the cell nucleus required to eradicate a    tumour cell population. Phys Med Biol. 2006; 51(9).-   R19. Huang C Y, Guatelli S, Oborn B M, Allen B J. Microdosimetry for    targeted alpha therapy of cancer. Comput Math Methods Med. 2012;    2012.-   R20. Miller B W, Frost S H L, Frayo S L, Kenoyer A L, Santos E,    Jones J C, et al. Quantitative single-particle digital    autoradiography with α-particle emitters for targeted radionuclide    therapy using the iQID camera. Med Phys. 2015; 42(7):4094-105.-   R21. Akabani G, Kennel S J, Zalutsky M R. Microdosimetric analysis    of α-particle-emitting targeted radiotherapeutics using histological    images. J Nucl Med. 2003; 44(5):792-805.-   R22. Humm J L, Roeske J C, Fisher D R, Chen G T Y. Microdosimetric    concepts in radioimmunotherapy. Med Phys. 1993; 20(2):535-41.-   R23. Hindorf C, Emfietzoglou D, Linden O, Kostarelos K, Strand S E.    Internal microdosimetry for single cells in radioimmunotherapy of    B-cell lymphoma. Cancer Biother Radiopharm. 2005; 20(2):224-30.-   R24. Elbast M, Saudo A, Franck D, Petitot F, Desbree A.    Microdosimetry of alpha particles for simple and 3d voxelised    geometries using MCNPX and Geant4 monte carlo codes. Radiat Prot    Dosimetry. 2012; 150(3):342-9.-   R25. Ma W, Wang X, Liu W, Ma H, Su Y, Yang Y, et al. A Theoretical    Model for Predicting and Optimizing In Vitro Screening of Potential    Targeted Alpha-Particle Therapy Drugs. Radiat Res. 2019; 191(5):475.-   R26. Chouin N, Bernardeau K, Davodeau F, Chérel M, Faivre-Chauvet A,    Bourgeois M, et al. Evidence of Extranuclear Cell Sensitivity to    Alpha-Particle Radiation Using a Microdosimetric Model. I.    Presentation and Validation of a Microdosimetric Model. Radiat Res.    2009; 171(6):657-63.-   R27. Akabani G, McLendon R E, Bigner D D, Zalutsky M R. Vascular    targeted endoradiotherapy of tumors using alpha-particle-emitting    compounds: Theoretical analysis. Int J Radiat Oncol Biol Phys. 2002;    54(4):1259-75.-   R28. Charlton D E. Radiation effects in spheroids of cells exposed    to alpha emitters. Int J Radiat Biol. 2000; 76(11):1555-64.-   R29. ICRU. Report 36. Microdosimetry. 1983.-   R30. Kellerer A M, Chmelevsky D. Concepts of microdosimetry—I.    Quantities. Radiat Environ Biophys. 1975; 12(2):61-9.-   R31. Kellerer A M. Analysis of Patterns of Energy Deposition. In:    Ebert H G, editor. Second Symposium on Microdosimetry. Stresa    (Italy); 1970. p. 107-36.-   R32. Kellerer A M, Rossi H H. A Generalized Formulation of Dual    Radiation Action. Radiat Res. 1978; 75(3):471-88.-   R33. Incerti S, Baldacchino G, Bernal M A, Capra R, Champion C,    Francis Z, et al. The Geant4-DNA project. Int J Model Simulation,    Sci Comput. 2010; 1(2):157.-   R34. Incerti S, Ivanchenko A, Karamitros M, Mantero A, Moretto P,    Tran H N, et al. Comparison of GEANT4 very low energy cross section    models with experimental data in water. Med Phys. 2010;    37(9):4692-708.-   R35. Bernal M A, Bordage M C, Brown J M C, Davídková M, Delage E, El    Bitar Z, et al. Track structure modeling in liquid water: A review    of the Geant4-DNA very low energy extension of the Geant4 Monte    Carlo simulation toolkit. Phys Medica [Internet]. 2015 December    [cited 2018 Sep. 14]; 31(8):861-74. Available from:    http://www.ncbi.nlm.nih.gov/pubmed/26653251-   R36. Incerti S, Kyriakou I, Bernal M A, Bordage M C, Francis Z,    Guatelli S, et al. Geant4-DNA example applications for track    structure simulations in liquid water: A report from the Geant4-DNA    Project. Med Phys. 2018; 45(8):e722-39.-   R37. Berger M J, Coursey J S, Zucker M A, Chang J. ESTAR, PSTAR and    ASTAR: Computer Programs for Calculating Stopping-Power and Range    Tables for Electrons, Protons, and Helium Ions (version 2.0.1)    [Internet]. NIST Standard Reference Database 124. 2017. Available    from: http://physics.nist.gov/Star-   R38. Guerard F, Gestin J F, Brechbiel M W. Production of    [211At]-astatinated radiopharmaceuticals and applications in    targeted a-particle therapy. Cancer Biother Radiopharm. 2013;    28(1):1-20.-   R39. Agostinelli S, Allison J, Amako K, Apostolakis J, Araujo H,    Arce P, et al. Geant4—a simulation toolkit. Nucl Instruments Methods    Phys Res Sect A Accel Spectrometers, Detect Assoc Equip [Internet].    2003 Jul. [cited 2018 Sep. 14]; 506(3):250-303. Available from:    http://linkinghub.elsevier.com/retrieve/pii/S0168900203013688-   R40. Allison J, Amako K, Apostolakis J, Araujo H, Arce Dubois P,    Asai M, et al. Geant4 developments and applications. IEEE Trans Nucl    Sci [Internet]. 2006 February [cited 2018 Sep. 14]; 53(1):270-8.    Available from: http://ieeexplore.ieee.org/document/1610988/-   R41. Allison J, Amako K, Apostolakis J, Arce P, Asai M, Aso T, et    al. Recent developments in Geant4. Nucl Instruments Methods Phys Res    Sect A Accel Spectrometers, Detect Assoc Equip [Internet]. 2016 Nov.    1 [cited 2018 Sep. 14]; 835:186-225. Available from:    https://www.sciencedirect.com/science/article/pii/S0168900216306957-   R42. Kellerer A M. Chord-Length Distributions and Related Quantities    for Spheroids. Radiat Res. 1984; 98(1):425-37.-   R43. Sgouros G, Roeske J C, McDevitt M R, Palm S, Allen B J, Fisher    D R, et al. MIRD pamphlet No. 22 (Abridged): Radiobiology and    dosimetry of α-particle emitters for targeted radionuclide therapy.    J Nucl Med. 2010; 51(2):311-28.-   R44. Paganetti H. Relative biological effectiveness (RBE) values for    proton beam therapy. Variations as a function of biological    endpoint, dose, and linear energy transfer. Phys Med Biol. 2014 Nov.    21; 59(22):R419-72.

1. A method comprising: determining volumetric data associated with anobject of interest, wherein the volumetric data comprise a plurality ofdomains; determining, for each of the plurality of domains, a pluralityof distributions for characterizing interactions of particles of aparticle beam with the corresponding domain, wherein the plurality ofdistributions comprises a first distribution indicative of energyimparted in a corresponding domain due to a particle traveling in thedomain; determining an analytical function based on one or more of theplurality of distributions; determining, based on the analyticalfunction and for each of the plurality of domains, first data comprisingenergy imparted by the particle beam to the corresponding domain; andoutputting, based on the first data, data associated with treatment ofthe object of interest by the particle beam.
 2. The method of claim 1,wherein the data associated with treatment comprises one or more of (1)a three-dimensional distribution of dose or (2) a segment-averaged,restricted, dose averaged linear energy transfer for the particle beam,and further comprising adjusting a treatment plan based on one or moreof the first data or the data associated with the treatment.
 3. Themethod of claim 1, wherein the first data comprises one or more of alinear energy transfer imparted by the particle beam to thecorresponding domain or a dose imparted by the particle beam to thecorresponding domain.
 4. The method of claim 1, wherein the particlebeam comprises a beam of one or more of protons, neutrons, positiveions, electrons, or alpha particles.
 5. The method of claim 1, whereindetermining the plurality of distributions comprises using one or moresimulations to generate the plurality of distributions.
 6. The method ofclaim 1, wherein the plurality of distributions comprises a seconddistribution indicative of segment length of a particle path in thedomain and a third distribution indicative of energy imparted in thedomain due to a collision of a particle, wherein the segment lengthcomprises a distance a particle of the particle beam is predicted totravel after entering the domain before coming to rest.
 7. (canceled) 8.The method of claim 1, further comprising determining, for a particleenergy of the particle beam and based on the plurality of distributions,an energetic kernel, wherein the energetic kernel comprises a firstaverage of at least one of the plurality of distributions and a firstvariance of the at least one of the plurality of distributions.
 9. Themethod of claim 8, further comprising performing a first convolution ofan energy fluence of the particle beam with the first average and aperforming a second convolution of the energy fluence of the particlebeam with the first variance, wherein the first data is determined basedon a result of the first convolution and the second convolution.
 10. Themethod of claim 1, wherein the data associated with treatment plan isbased on a model that accounts for one or more of variations of linearenergy transfer in a domain, variations of dose in a domain, variationsof segment length of paths of particles entering a domain, variations ofwhether particles come to rest in a domain, variations in a number ofcollisions of a particle in a domain, or variations in an amount ofenergy imparted in a collision of a particle in a domain.
 11. The methodof claim 1, wherein the volumetric data comprises one or more ofgeometric data associated with the object of interest, data comprising aplurality of voxels, data associated with a cell of the object ofinterest, data associated with a tissue of the object of interest, ordata associated with a macroscopic structure of the object of interest.12. The method of claim 1, wherein the volumetric data comprises one ormore of a geometrical model or a spatial distribution indicative of theobject of interest.
 13. The method of claim 1, wherein the volumetricdata is generated based on imaging data associated with the object ofinterest.
 14. The method of claim 1, wherein the domains comprisesubdivisions of a biological structure.
 15. The method of claim 1,wherein one or more of the plurality of domains vary in one or more ofshape, size, or arrangement to represent corresponding biologicalfeatures of the object of interest.
 16. (canceled)
 17. (canceled) 18.The method of claim 1, wherein determining the analytical functioncomprises determining the analytical function based on the firstdistribution.
 19. The method of claim 1, wherein determining theanalytical function comprises determining the analytical function basedon fitting the one or more of the plurality of distributions to theanalytical function.
 20. The method of claim 1, wherein the plurality ofdomains of the volumetric data indicate a shape of one or more of acell, a nucleus of a cell, or a tissue of the object of interest. 21.The method of claim 1, wherein the plurality of domains of thevolumetric data indicate an arrangement of one or more of a cell, anucleus of a cell, or a tissue within the object of interest. 22-44.(canceled)
 45. A device comprising: one or more processors; and memorystoring instructions that, when executed by the one or more processors,cause the device to: determine volumetric data associated with an objectof interest, wherein the volumetric data comprise a plurality ofdomains; determine, for each of the plurality of domains, a plurality ofdistributions for characterizing interactions of particles of a particlebeam with the corresponding domain, wherein the plurality ofdistributions comprises a first distribution indicative of energyimparted in a corresponding domain due to a particle traveling in thedomain; determine an analytical function based on one or more of theplurality of distributions; determine, based on the analytical functionand for each of the plurality of domains, first data comprising energyimparted by the particle beam to the corresponding domain; and output,based on the first data, data associated with treatment of the object ofinterest by the particle beam.
 46. A system comprising: a particle beamgenerator; and at least one processor communicatively coupled to theparticle beam generator and configured to: determine volumetric dataassociated with an object of interest, wherein the volumetric datacomprise a plurality of domains; determine, for each of the plurality ofdomains, a plurality of distributions for characterizing interactions ofparticles of a particle beam from the particle beam generator with thecorresponding domain, wherein the plurality of distributions comprises afirst distribution indicative of energy imparted in a correspondingdomain due to a particle traveling in the domain; determine ananalytical function based on one or more of the plurality ofdistributions; determine, based on the analytical function and for eachof the plurality of domains, first data comprising energy imparted bythe particle beam to the corresponding domain; and output, based on thefirst data, data associated with treatment of the object of interest bythe particle beam.